Consequences of paternally inherited effects on the genetic evaluation of maternal effects

Mixed models are commonly used for the estimation of variance components and genetic evaluation of livestock populations. Some evaluation models include two types of additive genetic effects, direct and maternal. Estimates of variance components obtained with models that account for maternal effects have been the subject of a long-standing controversy about strong negative estimates of the covariance between direct and maternal effects. Genomic imprinting is known to be in some cases statistically confounded with maternal effects. In this study, we analysed the consequences of ignoring paternally inherited effects on the partitioning of genetic variance. We showed that the existence of paternal parent-of-origin effects can bias the estimation of variance components when maternal effects are included in the evaluation model. Specifically, we demonstrated that adding a constraint on the genetic parameters of a maternal model resulted in correlations between relatives that were the same as those obtained with a model that fits only paternally inherited effects for most pairs of individuals, as in livestock pedigrees. The main consequence is an upward bias in the estimates of the direct and maternal additive genetic variances and a downward bias in the direct-maternal genetic covariance. This was confirmed by a simulation study that investigated five scenarios, with the trait affected by (1) only additive genetic effects, (2) only paternally inherited effects, (3) additive genetic and paternally inherited effects, (4) direct and maternal additive genetic effects and (5) direct and maternal additive genetic plus paternally inherited effects. For each scenario, the existence of a paternally inherited effect not accounted for by the estimation model resulted in a partitioning of the genetic variance according to the predicted pattern. In addition, a model comparison test confirmed that direct and maternal additive models and paternally inherited models provided an equivalent fit. Ignoring paternally inherited effects in the maternal models for genetic evaluation can lead to a specific pattern of bias in variance component estimates, which may account for the unexpectedly strong negative direct-maternal genetic correlations that are typically reported in the literature.


Background
Genetic evaluation based on mixed models [1] has become the main tool for prediction of breeding values in livestock populations. The most commonly used parameterization accounts for the direct polygenic additive genetic effect inherent to each individual and for several systematic effects (e.g., sex, herd, season, etc.), as well as a residual source of variation. However, some traits may also be affected by maternal genetic effects that refer to the influence from the genome of the dam (e.g., milk production of the dam). These maternal genetic effects are commonly considered in genetic evaluation models for some livestock populations [2,3]. Nevertheless, variance components estimates that are obtained with models that account for maternal effects have been the subject of a long-standing controversy about the unexpectedly strong negative estimates of the covariance between direct and maternal genetic effects [4][5][6].
Genomic imprinting [7,8] is the total or partial silencing of paternal or maternal alleles in the progeny. Several specific genes with imprinted inheritance have been identified in livestock species [9], such as the callipyge mutation in sheep [10,11] and the insulin-like growth factor 2 (IGF2) gene in pig [12,13]. Imumorin et al. [14] also reported that genomic imprinting is involved in the determinism of quantitative trait loci (QTL) in cattle. More generally, genomic imprinting has been described as a widespread phenomenon in other mammals [15]. Studies on pig growth [16] and bovine fatness traits [17] have also provided statistical evidence for the existence of genetic variation associated with paternal parent-oforigin inheritance that can account for over 10 % of the phenotypic variability in these traits. More recently, Neugenbauer et al. [18,19] and Tier and Meyer [20] confirmed this result on other cattle and pig traits.
Maternal effects can be statistically confounded with maternal inheritance, as noted by Hager et al. [21]. Furthermore, Meyer and Tier [22] suggested that ignoring parent-of-origin effects may be one of the possible causes for the unexpectedly large negative estimates of the correlation between direct and maternal genetic effects. In this paper, we provide further evidence that supports this latter argument by showing how the existence of a paternal parent-of-origin inheritance caused by the silencing of the alleles in the maternal gametes by imprinting can bias estimates of variance components in models that account for maternal effects. Specifically, we show that for any two individuals that were not related within the paternal lineage, direct and maternal additive genetic effects may generate the same expected covariance between relatives as paternal parent-of-origin or gametic effects. For most livestock pedigrees, the consequence of this is a downward bias in the estimate of the direct-maternal genetic covariance when paternal inheritance occurs and is not properly accounted for. To confirm this postulate, we developed a simulation experiment. In addition, we analysed a dataset of weight records at 210 days for the Pirenaica beef cattle population and we discuss the results in the light of this new evidence.

Quantitative genetic models
Consider a model with direct and maternal additive genetic effects. Hereafter, we use the subscripts "o" and "m" to differentiate between direct and maternal components, respectively. In addition, capital letters are used to identify different individuals. Under this model, the phenotypic values of individuals A and B are defined as: where μ is the general mean, a oX and a mX are the direct and maternal additive genetic effects, respectively, and e X is the residual for individual X = {A, B, A D , B D }. According to this notation, A D and B D stand for the dams of A and B, respectively.
Information for the estimation of variance components estimation originates from covariances between relatives. Under the model with maternal effects described above, and assuming independence of the residuals, the covariance between the phenotypic values of individuals A and B is: where σ a 2 , σ m 2 , and σ am are the direct additive genetic variance, the maternal additive genetic variance, and the direct-maternal genetic covariance, respectively, and r(X,Y) represents the additive relationship between individuals X and Y. In this context, the main sources of information to estimate maternal parameters are the observed correlations between phenotypes of full sibs (FS), maternal half sibs (MHS) and dam-offspring (DO) pairs. Equated to their expectations under the maternal model, the following expressions are obtained: Now, let: A little algebra shows that under condition (1) the previous expressions reduce to: This is exactly the expectation of the covariance between relatives for these same relationships under a paternal gametic model of inheritance. In the Appendix section, we show that this result can be generalised to any genealogical relationship between a pair of individuals, on the condition that one was not an ancestor of the other within the paternal lineage. Specifically, the covariance between the phenotypic values of any two individuals (say A and B) under the condition set by Equation (1)  where A S and B S are the sires of A and B, respectively. This striking result implies that the parameterization of the genetic covariance under the standard maternal animal model is confounded with a model that explains the whole genetic variation as paternally inherited; i.e. covðy A ; y B Þ ¼ r g ðA P ; B P Þσ 2 s ; where r g (A p ,B p ) is the gametic relationship [23] between paternal gametes of individuals A and B (A p and B p ) and σ s 2 is the variance caused by paternally inherited gametic effects. Stated the other way around, the existence of paternally inherited variation that is not accounted for by the maternal model will tend to be assigned to the remaining variance components according to Equation (1) in order to accommodate the observed correlations between relatives. As a consequence, variance component estimation will explicitly inflate estimates of the direct and maternal additive genetic variances, with a magnitude close to twice and one half the paternally inherited variance, respectively, and it will reduce the direct-maternal genetic covariance by a value equal to the paternally inherited variance. The following section describes a simulation experiment that we developed in order to test this result.

Simulation study
Each simulated dataset consisted of a base population of 1000 individuals (500 males and 500 females) and two generations of 5000 phenotyped individuals that comprised 2500 males and 2500 females. Phenotypic records were generated for the individuals of the last two generations based on the following model: where y is the vector of phenotypic data, μ is the general mean, set to 100 units, a o and a m are vectors of the direct and maternal genetic effects, s is a vector of paternally inherited gametic effects, and e is the vector of residuals, while Z a , Z m and Z s are incidence matrices. Covariances between genetic (a o and a m ) and gametic effects (s) were assumed to be zero and where I is the identity matrix and σ e 2 is the residual variance. The variances of the direct and maternal additive genetic effects are: V ar a o a m ¼ A⊗T; where A is the numerator relationship matrix and For the paternally inherited gametic effects (s): where G is the gametic relationship matrix Five simulation scenarios were developed: Scenario 1 corresponds to a pure direct additive model of inheritance with the following simulated parameters: Scenario 2 assumes that the only source of genetic variation is the paternally inherited variance: Scenario 3 combines both sources of genetic variation, direct additive and paternally inherited variance. The simulated parameters were: Scenario 4 corresponds to the covariance structure that is assumed in a standard maternal animal model. A strong negative covariance between direct and maternal additive genetic effects was also simulated in order to mimic estimates that are frequently achieved in livestock populations for maternally influenced traits: Finally, Scenario 5 includes all three types of genetic effects: direct and maternal additive and paternally inherited. The simulated parameters were: For each of these simulation scenarios, a total of 10 independent populations were generated.

Pirenaica beef cattle data
In addition to the simulated scenarios, we analysed data for phenotypic records on the Pirenaica beef cattle breed. The Pirenaica breed is a meat-type beef population from northern Spain with an approximate census of 20 000 individuals that are typically reared under extensive conditions [24]. The dataset consisted of 17 069 records for weight at 210 days (W210), with an average value of 266 kg and a raw standard deviation of 52.9 kg. Phenotypic records were adjusted according to the recommendations of the Beef Improvement Federation [25]. In addition, a pedigree file including 125 974 individual-sire-dam records was used. This information was provided by the National Breeders Confederation (Confederación Nacional de Asociaciones de Ganado Pirenaico-CONASPI; http://www.conaspi.es). Ethical approval for animal care and use was not required for this study since all data was field data from the Yield Recording System of the Pirenaica breed; furthermore, data was recorded by the stockbreeders themselves, under standard farm management, with no additional requirements.
An additional simulation study was performed based on the same scenarios as described in the previous section, but using genealogical information of the Pirenaica beef cattle breed and replacing each phenotypic record by a simulated value. The aim was to replicate the same structure of relationships between individuals of the real dataset.

Statistical models of estimation
Although a realistic implementation of all potential genetic sources of variation in the same hierarchical mixed linear model may result in non-estimable functions, the description of the full model is essential to understand the simplified parameterizations that are used below. The full statistical model is: where y is the vector of phenotypic data, b is the vector of systematic effects, consisting of the general mean for the simulated dataset and the two sexes and 16 age-groups of dam age of parity for the Pirenaica beef cattle data, and a o , a m , s and d are vectors of direct additive genetic, maternal additive genetic, paternal gametic and maternal gametic effects, respectively. Finally, p and h are vectors of random permanent maternal (9224 levels) and herd-year-season (2781 levels) environmental effects (not included in the simulated datasets), e is a vector of residuals, and X, Z a , Z m , Z s , Z d , Z p , and Z h are incidence matrices that link the effects with the phenotypic data. Variances and covariances for random sources of variation were defined as follows: where σ p 2 and σ h 2 represent the variances of the permanent maternal and herd-year-season effects, σ d 2 is the variance of maternal gametic effects, and σ sd is the covariance between paternal and maternal gametic effects. All other parameters were defined previously.
Note that under a pure direct additive model: A number of models based on reduced parameterizations of this full model were fitted to both the simulated and the real data: Model A: y = Xb + Z a a o + e, Model S: y = Xb + Z s s + e, Model AS: y = Xb + Z a a o + Z s s + e, Model SD: y = Xb + Z s s + Z d d + e, Model AM: y = Xb + Z a a o + Z m a m + e, Model AMS: y = Xb + Z a a o + Z m a m + Z s s + e, When the Pirenaica beef data was used, all models also included the terms z p p and z h h.
For statistical analysis, we applied a Bayesian approach through a Gibbs sampler [26]. The prior distributions for systematic effects and variance components were assumed to be bounded uniform, and the prior distributions for the genetic effects were multivariate Gaussian distributions with mean zero and variance as defined in the previous section. Analyses were run for each simulated or real dataset as a single chain of 525 000 cycles with the first 25 000 iterations being discarded. Convergence was checked by visual inspection of the chains and by applying the test of Raftery and Lewis [27]. All samples were stored to calculate summary statistics.

Model comparison
Models were compared using the pseudo log-marginal probability of the phenotypic data [28]. If we consider the data vector y = (y i, y -i ), where y i is the i th datum and y -i is the vector of data with i th datum deleted, the conditional predictive distribution has a probability density equal to: where θ is the vector of parameters. Therefore, p(y i |y -i ) can be interpreted as the probability of each datum given the rest of the data and is known as the conditional predictive ordinate (CPO) for the i th datum. The pseudo log-marginal probability of the data (LogCPO) is then: The collection of conditional predictive densities is equivalent to the marginal probability of the data [29]. A Monte Carlo approximation of the LogCPO suggested by Gelfand [28] is: and N is the number of Markov chain Monte Carlo (McMC) draws, and θ j is the j th draw from the posterior distribution of the corresponding parameter. The interpretation of the results of this test is that the larger is LogCPO, the higher is the marginal probability of data and, thus, the better is the relative fit.

Scenario 1: additive model of inheritance
In Scenario 1, all genetic variation was generated as direct additive; there were no maternal genetic or paternally inherited effects. The average posterior mean estimates of variance components and ratios of components of variance are in Table 1. When Model A was used, the estimate of the direct additive genetic variance almost matched the simulated value. Furthermore, when the estimation model only included paternal gametic effects (Model S), a significant proportion of the variation was assigned to these effects and the remaining additive variation was incorporated into the residual variance. Implementation of a more complex parameterization, such as Model AS, attributed a small, although negligible, amount of the simulated variability to the paternally inherited effects. As expected, Model SD absorbed close to half of the additive genetic variance of the paternal gametic effects and another half of the maternal effects. In addition, the correlation between paternal and maternal effects was close to 1. Finally, the models that included an additive maternal genetic effect (AM and AMS) also correctly assigned most genetic variation to the direct additive effects, and only small amounts to other variance components. The results of the model comparison test showed that Models A and SD presented the best fit, followed closely by Models AM, AS and AMS. Model S had the worst fit.

Scenario 2: paternal inheritance model
Scenario 2 (Table 2) only included a variance of paternally inherited effects of 250, without considering any other source of genetic variance. The average posterior mean estimate of variance due to paternally inherited effects under Model S was close to the simulated value. However, with Model A, estimates of direct additive genetic variance were on average close to 150, although this variance component was not included in the simulation process. Model AS assigned most of the variance to paternally inherited effects and a small amount to the direct additive effect. Model SD assigned most of the genetic variance to paternal gametic effects, although negligible amounts of variability were also assigned to the gametic maternal variance component. Model AM, which corresponds to the standard maternal animal model, strongly overestimated the direct and maternal additive genetic variances, and resulted in a negative estimate of the covariance between direct and maternal additive genetic effects, which resulted in a negative correlation close to −1. This tendency was reduced when paternally inherited effects were included in the model (Model AMS). Average posterior mean estimates of variances of both direct and maternal additive genetic effects and the covariance between them were smaller than in Model AM. However, the genetic correlation was still negative (−0.90) and a substantial percentage of variability was assigned to paternally inherited effects. It should be noted that the residual σ a 2 and σ m 2 are the direct and maternal additive genetic variances and σ am is the direct-maternal genetic covariance. σ s 2 and σ d 2 are the paternal and maternal gametic variances and σ sd is the covariance between the paternal and maternal gametic effects. Finally, σ e 2 is the residual variance and L is the average logCPO across 10 replicates.
variance was inflated in Model A and smaller than the simulated value in Model AM. The model with the best fit was Model S, which mimicked the simulation model, and was closely followed by Models AS and SD, which included paternally inherited effects, and by Model AM. Model A had the worst fit.  Table 3. The first and second models of estimation (Models A and S) clearly overestimated the additive genetic variance and the paternal gametic variance, respectively. This overestimation disappeared with Model AS, which replicated the model of simulation by including a paternal gametic effect, and the estimates were close to the simulated values. Models AM and AMS overestimated the additive genetic variance and assigned some of the variance to the maternal genetic effects. Furthermore, these models generated a negative covariance between additive and maternal genetic effects. This overestimation of variances was reduced with Model AMS, for which a substantial amount of the variability was assigned to the paternally inherited effects, although less than the simulated value. The model comparison test indicated that Models AS, SD, AM and AMS had the best fit, with similar average LogCPO values, whereas Models A and S had a substantially poorer fit.

Scenario 4: maternal animal model
Scenario 4 corresponds to the standard maternal animal model. Direct and maternal additive genetic variances of 500 and 250 were generated, but with a strong negative correlation between these two effects (−0.707). Average posterior mean estimates are in Table 4. As expected, the average posterior mean estimates from Model AM σ a 2 and σ m 2 are the direct and maternal additive genetic variances and σ am is the direct-maternal genetic covariance. σ s 2 and σ d 2 are the paternal and maternal gametic variances and σ sd is the covariance between the paternal and maternal gametic effects. Finally, σ e 2 is the residual variance and L is the average logCPO across 10 replicates. σ a 2 and σ m 2 are the direct and maternal additive genetic variances and σ am is the direct-maternal genetic covariance. σ s 2 and σ d 2 are the paternal and maternal gametic variances and σ sd is the covariance between the paternal and maternal gametic effects. Finally, σ e 2 is the residual variance and L is the average logCPO across10 replicates.
almost matched the simulated values for the direct and maternal additive genetic effects, and the covariance between them. Moreover, all models except AM gave biased estimates of variance components. In particular, Model A resulted in an average posterior mean estimate of the direct additive genetic variance that was smaller than the simulated value. Model S assigned a significant proportion of the variation to paternally inherited effects, although they were not included in the simulation model. This was confirmed with Model AS, which assigned a significant proportion of the variation to the paternally inherited effects. Model SD assigned a significant proportion of the variance to both paternal and maternal gametic effects, but the covariance between them was almost null. Finally, with Model AMS, the presence of a significant amount of variation in the paternal gametic effects reduced the magnitude of the direct and maternal additive genetic variance components compared to the results obtained with Model AM. For this scenario, the models with the best fit were SD and AM, followed by Models AS and S. The models that fitted the simulated data worst were A and AMS.

Scenario 5: direct, maternal and paternal inheritance model
Scenario 5 included direct and maternal additive genetic, and paternal gametic effects. Direct and maternal additive genetic effects with variances of 500 and 250 were generated, with a strong negative correlation between them (−0.707). The paternally inherited variance was 250. The average posterior mean estimates are in Table 5. The average results of Model AMS were very close to the simulated values. However, the standard deviation of the estimates was substantially higher than for the other models and also than for the previous scenarios. The results of Model A were very close to the simulated value for the direct additive genetic effect (510.8 vs 500), although residual variance was significantly overestimated. Model S, AS and SD overestimated the paternal gametic variance, which captured the sources of variation that were ignored in the estimation model. It is remarkable that Model SD provided almost a null covariance between sire and dam gametic effects. Finally, the results of Model AM reflected a redistribution of the ignored paternally inherited variance by increasing the direct and maternal additive genetic variances (to 924.6 and 371.0, respectively), and by generating a strong negative directmaternal genetic covariance (−473.5). The results of the LogCPO test indicated that Model AMS, which replicated the simulation model, had the best fit. It was closely followed by Models AM and SD, whereas Models AS, S and A fitted the simulated data worst.

Pirenaica beef cattle breed data
The results of the analysis using real data on the Pirenaica beef cattle breed are in Table 6. With Model A, the posterior mean estimate of the direct heritability (h 2 d ) was equal to 0.367, but with Model AS, it was only equal to 0.174, whereas the posterior mean estimate of paternal gametic heritability was 0.379, which increased to 0.505 with Model S that did not take direct additive genetic effects into account. Model SD returned the highest posterior mean estimate of the paternal heritability (0.463), a lower value for maternal heritability (0.130) and a very low estimate of the direct-maternal genetic correlation (0.059). The standard maternal animal model (AM) produced a very large posterior estimate of direct additive heritability (0.471), along with a smaller posterior estimate of maternal heritability (0.157), and the correlation between direct and maternal effects was strong and negative (−0.734). Finally, with Model AMS, both heritabilities fell to very low values Table 4 Averages of posterior means (and posterior standard deviations) of (co)variance components for Scenario 4 (σ a 2 = 500, (0.045 and 0.056) and a significant proportion of the genetic variation was assigned to the paternal gametic effects (0.436). The results of the comparison of models based on LogCPO showed that Model SD had the best predictive ability, followed by Models AM and ASs. Model AMS was next, and Models S and A showed the worst fit.

Discussion
We present a theoretical development that shows that the structure of the (co) variance components assumed by the standard maternal animal model is partially confounded with a model that considers a paternal gametic effect. Meyer and Tier [22] postulated that this may be one  σ a 2 and σ m 2 are the direct and maternal additive genetic variances and σ am is the direct-maternal genetic covariance. σ s 2 and σ d 2 are the paternal and maternal gametic variances and σ sd is the covariance between the paternal and maternal gametic effects. σ h 2 and σ p 2 are the herd and maternal permanent environmental variances. σ e 2 is the residual variance. h a 2 , h m 2 , h s 2 and h d 2 are the ratios of the direct additive genetic, maternal additive genetic, paternal gametic and maternal gametic variances to the phenotypic variance, respectively. r am is the genetic correlation between direct and maternal additive effects and r sd is the correlation between paternal and maternal gametic effects. c h 2 and c p 2 are the ratios of herd and permanent maternal environmental effects to the phenotypic variance. σ a 2 and σ m 2 are the direct and maternal additive genetic variances and σ am is the direct-maternal genetic covariance. σ s 2 and σ d 2 are the paternal and maternal gametic variances and σ sd is the covariance between the paternal and maternal gametic effects. Finally, σ e 2 is the residual variance and L is the average logCPO across 10 replicates possible cause for the strong negative estimates of the genetic correlation between direct and maternal genetic effects that are typically obtained. We derived results that show that, for most of the genetic relationships in livestock pedigrees, the covariance between relatives generated by a paternally inherited effect is mimicked by the maternal animal model when variance components follow the pattern defined by Equation (1). Thus, when paternally inherited effects are ignored in the estimation model, its variance is assigned to the remaining variance components in order to accommodate the observed correlations between relatives. Under the maternal animal model, its consequences are an inflation of the direct and maternal additive genetic variances, and a downward bias of the direct-maternal genetic correlation.
These results were confirmed by the simulation study. When simulation and estimation models were analogous, the Bayesian analysis was able to provide suitable estimates of the variance components. Furthermore, when the true simulation structure of the (co)variance components was included in the estimation model together with other variance components, the procedure correctly estimated the true parameters, whereas estimates of the non-simulated variance components were small. However, the presence of paternal gametic effects (Scenarios 2 and 3) and their absence in the estimation model resulted in an overestimation of the direct additive genetic variance, as previously described by other authors [17,22]. Furthermore, when paternally inherited effects that play a significant role in the genetic variation were ignored in a model that included correlated direct and maternal additive genetic effects, the direct additive genetic variance component was increased by almost twice the paternal gametic variance, whereas the genetic covariance between direct and maternal effects became negative and the estimate of the maternal additive genetic variance increased by almost one half of the simulated paternal gametic variance, as was observed in Scenarios 2, 3 and 5. As described above, this redistribution of the gametic paternal variance is linked to the constraint imposed in Equation (1). These results were also consistent with the presence of a significant percentage of paternal gametic variance for Models S, AS, SD and AMS in Scenario 4, where only direct and maternal additive effects were simulated, and with the increase of paternal gametic variance in Models S, AS and SD in Scenario 5. Finally, the decrease in covariance between paternal and maternal effects under Model SD in Scenarios 4 and 5 may also be attributed to this phenomenon.
More generally, when any relevant variance component in the simulation model was not considered in the estimation, substantial biases on the variance component estimation were noted. For example, ignoring a true direct additive genetic effect while assuming paternal and maternal gametic effects in the estimation model, led to the absorption of close to half of the additive genetic variance by the paternal gametic variance and half by the maternal gametic variance. In addition, the increase in the estimate of the covariance between paternal and maternal gametic effects was also of this magnitude. This distribution of the direct additive variance followed the pattern described in Equation (2) and was exploited in the model proposed by Neugenbauer et al. [18,19] in order to detect paternal and maternal inherited effects.
The model comparison test that we performed provided support to our main argument. To illustrate this point, in Scenarios 2 and 3, where paternally inherited effects were included in the simulation model, the goodness of fit of Model AM was very close to that obtained with Models S and AS, which mimicked the respective simulation models. Moreover, under Scenarios 4 and 5, Model SD performed almost as well as Models AM or AMS, which replicated the true simulation models. Clearly, the constraint imposed by Equation (1) explains the similarity observed between the maternal animal model and the paternally inherited model in those cases.
It should be emphasized that the constraint imposed by Equation (1) is not a pure parametric equivalence, since it is limited by the condition that available correlations between relatives do not involve genealogical relationships where one individual is an ancestor of the other within the paternal lineage. Therefore, the ability to separate direct and maternal additive genetic effects and sire imprinted effects can only be achieved if relevant records from sires and their offspring are available or, more generally, if there exist records from individuals that are genetically linked through the paternal pathway. However, in beef cattle datasets, phenotypes on sire-offspring pairs are usually rare. In the Pirenaica beef cattle dataset analyzed here, there were 527 463 paternal half-sibs, 15 128 maternal half-sibs, 2707 full-sibs, 2411 dam-progeny and only 4566 sire-progeny pairs with phenotypic information on both individuals. Thus, most of the information available for estimation of variance components was provided by half-sib relatives, whose covariance generated by paternally inherited effects can be also achieved by the maternal animal model under the constraint imposed by Equation (1). This latter explanation helps to understand the results of the analysis of the Pirenaica beef cattle dataset. Estimates of variance components obtained from Model A were concordant with previous results obtained in the same [30] and other beef cattle populations [31,32]. In turn, alternative estimation models, such as S, AS and SD, suggested that paternal gametic effects may play an important role in the variation of weight at 210 days in this breed. In contrast, under Model AM, maternal additive effects were also very relevant. Furthermore, when both effects were fitted in the estimation model (AMS), the maternal effects contributed less to the total variation than the paternal gametic effects. A similar reduction in the maternal contribution on weaning weight was also found in the Bruna dels Pirineus beef cattle breed when paternal gametic effects were considered in the estimation model [33]. Indeed, more evidence for the presence of parent-of-origin variation is derived from the fact that when Model SD was used, a very low correlation between paternal and maternal gametic effects was obtained. This is in contrast to the results reported by Neugenbauer et al. [18,19], which provided a high correlation between them, as in the outcome of the first simulation scenario, where a pure additive model was simulated. However, as observed with Scenarios 4 and 5, independence between paternal and maternal gametic effects, together with the presence of a significant paternally inherited variance component can also be due to the existence of a strong and negative genetic correlation between direct and maternal genetic effects. Given the structure of the information in the Pirenaica dataset, the genetic variance is distributed between paternal gametic and direct and maternal effects depending on the model of estimation and, thus, it is very difficult to clearly identify the true sources of variation. However, at this point, it is worth mentioning that the results of the simulation study using the Pirenaica's breed genealogical structure were similar to those obtained with the above described simulation study (results in Tables S1, S2, S3, S4 and S5 [See Additional file 1: Tables S1-S5]). This allowed us to rule out the possibility of a statistical artefact caused by the structure of the information. Nevertheless, the model comparison using logCPO indicated that Model SD was most suited to the Pirenaica beef cattle data. It should be noted that this model performed well in all simulation scenarios (Tables 1, 2, 3, 4 and 5), since it can accommodate a significant amount of paternal inheritance and also a strong negative correlation between direct or maternal effects, as in Scenario 4. In fact, Model AM had the second best predictive ability based on the LogCPO approach and showed a better fit than Model AS. However, the very large difference in LogCPO between Models SD and AM (−86178 vs −86378) and the results of Model AMS suggest that the existence of significant paternal gametic effects is plausible in the Pirenaica population.
Many studies have reported negative estimates of the genetic covariance between direct and maternal additive effects for beef cattle populations [4][5][6]. Several explanations have been proposed, including environmental correlations between progeny performance and maternal environmental effects [34,35] and sire × herd interactions [36,37]. An additional explanation that should be considered is the presence of paternally inherited effects, as previously suggested by Meyer and Tier [22]. Our evidence is consistent with the effect associated with imprinted genes, such as the insulin-like growth factor 2 (IGF2) gene involved in fetal growth [38,39] or postnatal growth and fat deposition in pigs [11,12,40] and mice [41]. Genomic imprinting of the IGF2 gene has also been shown to affect meat quality traits in beef cattle [42,43] and some mutations in IGF2 have been associated with carcass traits [44] and body weight [45].
Furthermore, a wealth of evidence on the role of parental inheritance in fetal and early growth in mammals supports the imprinting hypothesis [46,47]. This is related with parent-offspring conflict theories about the evolutionary origin of paternal effects [48], which are based on the hypothesis that the self-interest of the paternal genome is to increase the growth of the fetus or young individuals, whereas the aim of the maternal genome is to sufficiently limit growth to ensure the success of future offspring without jeopardizing the mother's survival. According to this hypothesis, genes that enhance growth tend to be more expressed if they are paternally inherited than when they are maternally inherited. This phenomenon should be even more marked in species with a low rate or absence of full siblings and a high frequency of multiple paternity [49], as is the case for cattle.
Finally, other causes that can explain overestimation of paternal gametic variance are a large difference in allele frequencies between sires and dams or the presence of significant dominance variation. However, these are not relevant for close populations and that include few fullsib relationships, respectively. In future studies, a very promising alternative to distinguish between paternally and maternally inherited effects arises from the use of molecular information that may contribute to better discriminate between paternal and maternal alleles.

Conclusions
Ignoring paternally inherited effects in genetic evaluations with maternal animal models can lead to a specific pattern of bias in the estimates of variance components. This phenomenon may account for the unexpectedly strong negative direct-maternal genetic correlations that are reported in the literature. In addition, we showed that the existence of paternal effects may play a role in the explanation of the partitioning of genetic variation in growth traits and possibly in some other relevant traits in livestock production. The main implication of this study is that the presence of relevant paternally inherited variation should be investigated when the maternal animal model provided negative correlation estimates between direct and maternal effects. In addition, in those cases, applying maternal animal models should be done with caution.