The efficiency of multi-generation selection on maternal traits, with implications for reindeer

Maternally affected traits, such as juvenile growth and survival, provide resilience in mammal species, in particular for reindeer living in extreme northern habitat. The genetic variation in such traits is caused by direct and maternal genetic effects (DGE and MGE, respectively). We used Willham’s variance-component approach and extended a family index with the focal individual and its fulland half-sibs to an approximated BLUP (pseudo-BLUP) by including the parents’ estimated breeding values. Most of the deviations of the predicted responses from the simulated ones were 4.1% for DGE and 5.3% for MGE. The benefits of index and BLUP selection are high in the case of negative correlation, large full-sib family and in particular, when maternal half-sibs are available. Higher economic value for MGE than for DGE is needed, since with equal heritabilities and economic weights for the effects the maternal response is 40 to 70% of the direct one. With negative correlation, records on collateral relatives beyond sibs are possibly needed. They would support also the prediction of MGE in uniparous reindeer lacking full-sib information.


Introduction
Maternal effects are common and important in all mammals and describe the effects a mother (dam) has on her offspring's phenotype. In general, such effects are social interaction or indirect genetic effects with possible environmental component (Dickerson 1947, Willham 1963, Moore et al. 1997, Bijma 2006, Wolf and Wade 2009. In natural populations maternal care and other maternal effects are important in determining the survival and other fitness traits of the young (Donohue 1999, Hereford and Moriuchi 2005, Räsänen and Kruuk 2007. In heterogeneous and extreme northern environments utilized by semi-domesticated reindeer with the habitat of the animal undergoing yearly changes, maternal effects have an important role in phenotypic plasticity providing a fast, adaptive (but temporary) response to the changing environment shared by offspring and mother (Mousseau and Fox 1998).
Animal breeding research has been occupied by the genetic basis of maternal effect as an indirect genetic effect (IGE) affecting the total genetic variation available for selection (Bijma et al. 2007). In meat producing livestock species, maternal effects are important for traits relating to growth such as birth weight, weaning weight, average daily gain etc.
Many maternal effects consist of behavior traits and one option would be to perform selection on behavioral observations (Chiang et al. 2002). In circumventing laborious data recording problems, animal breeders have facilitated the selection on maternal traits by the use of statistical methods. The variation in maternal traits can be modelled by Willham's (1963Willham's ( , 1972 variance component model expressing the phenotype as the sum of individual's direct effect and its' dam's maternal effects. Our goal is the development of methodology to assess the design and volume of data collection to reach satisfactory accuracy of breeding values for direct and maternal effects with multi-generation information. A general (analytical) method is desirable because the alternative, stochastic simulation, would give answers only to specific situations. Simulation have been used for example by Lourenco et al. (2013) to study the effect of genetic evaluation method and Maiorano et al. (2019) to study selection when using pooled semen. Where available, maternal half-sibs are expected to be efficient in predicting MGE (Willham 1963, Noble et al. 2014, Maiorano et al. 2019). Wray and Hill (1989) introduced an analytical method for approximating the use of multi-generation pedigree information for predicting the response to selection (often now called pseudo-BLUP). The method starts from a selection index on collateral relatives' phenotypic information and mimics the influence of pedigree information using parents' estimated breeding values (EBV's) and iterating over generations. The compact method lends itself to extensions, e.g. Villanueva et al. (1993) extended the method to multi-trait context, which was applied by Mulder and Bijma (2005) to study genetic evaluation in traits under GxE interaction. Dekkers (2007) used an Manuscript received September 2020 approximation to study the benefits of computing and using genomic EBV's. Mulder (2016) resorted to pseudo-BLUP approach in investigating the gains due to genomic selection for selection strategies under GxE interaction.
The objectives of our study are to: Extend and validate a pseudo-BLUP approximation to traits influenced by maternal effects. Study the outcome of selection schemes for direct and maternal effects with different economic weights, design and volume of data collection. We discuss implications for data collection and selection schemes, in particular for reindeer.

Materials and methods
Quantitative genetic framework A maternally affected phenotype ( ) is made of an individual's direct genetic effect (DGE) and its dam's (throughout the text subscript d) maternal genetic effects (MGE) (covering both pre-and postnatal effects), and the trait is expressed as ( deviation from the population mean): (1) where is (possibly correlated) environmental components of both indirect and direct effects and c stands for the effect of common environment shared by the full-sibs. We assume an infinitesimal model with normally distributed genetic effects (Fisher 1918, Bulmer 1980). We use the notation and for direct and maternal genetic variance, respectively and for their covariance. Further is the variance due to common environment and is the variance due to other environmental effects.
The genetic (co)variances among the selected parents are affected by (possibly different) selection in the two sexes and the covariance between the effect and selection criterion (phenotype) ( and covariance between and DGE and MGE, respectively) (Pearson 1903, Bulmer 1971, Robertson 1977. For example, the covariance amongst selected dams is where , is the intensity of selection and is the truncation point on the standardized scale.  Selection does not affect the (co)variances , and (0 referring to the base generation) due to segregation. Therefore in the generation t the aforementioned covariance is (the bar refers to the mean across sexes). The response to selection (or the visible change in the phenotype) after a generation of phenotypic selection (cf Willham 1972) is . An asymptotic response is reached within few generations when the Bulmer effect reaches a balance and , and . The response can be written , where and are the (co)variances among selected dams.

Selection index framework
We construct the family index based on the phenotype of an individual ( ), the mean phenotype of its full-sibs (excluding the individual) and the mean phenotype of its half-sibs (excluding the individual and the full-sib family where it belongs to).
of the individual is for either effect where b j is the respective index weight of the relative's information.
In a hierarchical mating structure, each sire (subscript s throughout the text) mates with f dams and each dam has n offspring. The variances of and in generation t are according to Willham (1963): and .
Usually a two-trait approach is used for maternal traits, while in constructing the selection index resorting to the same phenotypic information, we use the approach similar to Van Vleck (1970). The aggregate selection goal of maternally affected trait can be written H = [a m]w, where w is a 2 × 1 vector for economic weights for DGE and MGE, respectively (e.g. in the case of equal weighting, w is a vector of ones). The selection index is I = b x, where x is a vector of information sources (in addition to p, and also the dam's and sire's index and the mean index of half-sibs' dams) and b is a 6 × 1 vector of selection index coefficients solved from the ordinary selection index equation b = P -1 G w. The formulae for the covariances in the P matrix are: , , , and and the variance of index values amongst selected dams (and sires) is: and ) and dams of the half sibs is and and . The matrix P (in generation t) is a (co)variance matrix (2) of information sources without common environmental variance. (2) The G t matrix (for generation t) has a column for both direct and maternal effect corresponding to the covariance of a respective genetic effect with the sources of information in the index Alternatively we can construct the family index using independent information sources (within full-sib family deviations and full sibs within half-sibs, i.e. no covariances) and have , where the three independent information sources demonstrate the influence of full-sib information versus within family deviations and half-sib mean. The diagonal elements in P are: , and . The G matrix is: (4) Using independent information sources in the case of maternal half-sibs (S = sires per dam) the diagonal elements of P t are then: , and (with the off-diagonals being zero) and (5) Pseudo-BLUP The family index can be extended to include information available on the parents, i.e. their EBV, and the index for individual would be where is the mean of the paternal half-sibs' dams' EBV´s. We can mimic the effect of accrued pedigree information by repeating the calculations over generations and arrive at an approximated BLUP (or pseudo-BLUP) of breeding values (Wray and Hill 1979). In the pseudo-BLUP the available information for the index accumulates over generations, and the individual's EBVs become more accurate (corresponding to the effect of large and information rich relationship matrix in BLUP), which leads to higher response to selection. After the initiation of (co)variances, the steps for predicting the response at each generation for index selection and pseudo-BLUP are (i) calculation of index coefficients and responses, and (ii) updating (co)variances in P and G.
We used pseudo-BLUP method to study the influence of variable family full-sib family size and mating ratio. With both paternal and maternal half-sibs, different heritabilities and genetic correlation between DGE and MGE. In addition, the common environmental variance amongst full-sibs and variable economic weights for the effects was investigated.

Simulation
We compared the pseudo-BLUP prediction with stochastic simulation results using R-program (R Core Team 2019). The packages cpgen (Hauer 2015) and Matrix (Bates and Maechler 2018) were used for optimization of matrix calculations. Each new discrete generation was generated with 50 sires each mated with 10 females producing 10 offspring each with the selection intensity (proportion selected) across offspring generation being approximately 2.4 (2%) for males and 1.4 (20%) for females. High full-sib family size was chosen because selection for MGE is then more successful and the substantial differences in selection intensities in the sexes were used to mimic the possibilities in many species. At the beginning, we generated unrelated parents and grandparents, and by random mating, an offspring generation in order to produce maternally affected phenotypes and genetic values for the first selection candidates. We used the R-package MCMCglmm (Hadfield 2010) to generate the starting values for genetic effects in the founder population. Selected individuals were randomly mated to generate the offspring with consideration of the effect of possible inbreeding on the segregation variance utilizing the R-package pedigree (Coster 2012) in estimating the inbreeding coefficients. In the simulation the evaluation was based on either phenotype, family index or BLUP with full pedigree information. The selection was carried out for 10 generations with equal economic weights for the two effects, and an asymptotic response was reached. We ran all simulated cases multiple (50) times to establish reliable estimates (expressed with standard errors) for the means across replicates.
We tested the three cases with different genetic correlations ( ) between DGE and MGE (0.5, 0.0, -0.5) keeping the heritabilities for the effects at 0.30. We studied a whole range of genetic parameter values (results not shown) and here a high heritability for MGE was used to better demonstrate the potential selection response in MGE. The deviations due to the covariance of common environmental effects were not simulated.

Phenotypic, family index and pseudo-BLUP selection
Overall, the analytical methods predict in a satisfactory way the outcome of the simulated response to selection ( Fig. 1 and Table 2). The difference between the analytical prediction and the simulated response to selection was on average below 4.1% for DGE and 5.3% for MGE. Standard error of the simulated response mean in the last generation ranged from 0.003 to 0.007. For MGE, the predicted response to family index selection was slightly below the simulated gain for and vice versa for .
The family index and pseudo-BLUP are very accurate in predicting DGE. When , the response to phenotypic selection for MGE with = 0 is about half of that in DGE, reflecting the one generation gap in the expression of the maternal effect. With family index or BLUP, the difference is less drastic, but the response of DGE is still ∼70% higher than that of MGE. The phenotypic selection is the least efficient and does not yield any gain for < 0, while the family index and BLUP are qualitatively different with the latter yielding the highest response. With < 0, the MGE response to family index and BLUP selection is positive in contrast to the phenotypic selection and almost 40% of that in DGE.  The asymptotic rate of response is, in general, 75-80% of the first generation response but this difference is consistent over the cases, indicating that the latter could be used for comparing selection strategies. This has been observed in the earlier studies (cf Wray and Hill 1989).
With < 0 the pseudo-BLUP has a considerable prediction error (some 20%) indicating need to add information sources beyond full-and half-sibs and their parents' EBV's.

Variable family size and mating ratio
In general, with paternal half-sibs, higher f improves the response of DGE (with slight relative decrease for the maternal effect) while larger n improves the response of MGE (cf Table 4). With uncorrelated DGE and MGE, half-sib information alone does not capture any maternal effect variation, as implied by the G matrix in (4). Similar results have been found analytically by Lande and Kirkpatrick (1990), who did not consider the Bulmer effect or inbreeding.

Maternal half-sibs
While paternal half-sibs provide no information on MGE when = 0, maternal half-sibs are very valuable even then (Table 3) and with equal heritabilities the MGE response is twice the DGE one. Further when ≠ 0, the individual's own record is better in predicting MGE than its full-sibs.

Different heritabilities
When and is highly negative, then with paternal half-sib structure the negative response in MGE is moderated having data from large n ( Table 4, upper section). When is moderately negative, the MGE response is positive with large n. In general, when f increases the response for MGE effect decreases, if is small and positive (for example < 0.75, Table 4).
When (Table 4, lower section), and is highly negative (-0.75), selection favors MGE almost exclusively, even if f is high. When is less negative (-0.25) or positive, the response for DGE is higher as f increases and is, in general, lower for higher n. The selection response for MGE is again proportional to n and lower with larger f as was noticed for the equal heritability case. Table 2. The predicted (pred) and simulated (sim) rate of response (in σ p units of the base generation) to selection on phenotype, family index and approximated BLUP of the aggregate genotype of equally weighted direct ( ) and maternal ( ) effect. The response is shown for first generation and an average rate across generations 7 to 10. = 0.3 and = -0.5, 0.0, +0.5. The population is made of 50 sires each mated with 10 dams producing 10 offspring each and the selection intensity is 2.4 in males and 1.4 in females.

Economic weights
Equal weighting of DGE and MGE gives the highest total response independent of the sign of (Table 5A). When < 0 the response with prioritized MGE (w = [1 2]) gives a better MGE response, in particular for high . With > 0 , the response in MGE is substantial and close to that in DGE with moderate or high (results not shown). When either effect is ignored in the selection goal, the co-response in MGE smaller than the DGE one.

Common environmental variance
A high common environmental variance decreases only the accuracy of EBV for MGE (Table 5B). This is reflected also as a reduced overall response. Table 3. The effect of independent family information sources on selection response after one generation of selection when each dam is mated to multiple sires (maternal half sibs). f = 5 or 10 and n = 2 or 10, = 0.3 and = -0.5, 0.0, +0.5. Selection intensity was 2.06 in both sexes.

Methodology and results
The main result of the study is the development of an analytical approximation for an animal model BLUP selection for traits influenced by maternal effects. Reduction of genetic variance due to Bulmer effect is considered enabling the comparison for simulated truncation selection. Because the number of sires was 50 (implying negligible rate of inbreeding), in the pseudo-BLUP there was no need to adjust the segregation genetic variance for inbreeding, neither to reduce selection intensity due to finite population size (cf Meuwissen 1991). Khaw et al. (2014) mentioned the construction of pseudo-BLUP for IGE, while due to the large number of required parameters and subsequent complicated structures ultimately preferred to use simulation in their study. Maternal traits represent a simple case for IGE and in our approach only few additional parameters are needed in a single trait case to arrive at compact pseudo-BLUP expressions.
In general, BLUP selection gave (slightly) higher response compared to family index. The rate of genetic change asymptotes when selection is carried over several generations. The increased FS family size always benefits the selection of MGE. Although rarely available within the same generation, maternal half-sibs are naturally more desirable than paternal half-sibs (Maiorano et al. 2019) in increasing the accuracy of MGE predictions. In animal breeding studies is often found to be negative (Bradford 1972, Koch 1972, Robinson 1972, Bijma 2006. Phenotypic selection fails to achieve any gain for negative . In such a case, the family index and BLUP are qualitatively different. The cases with highly negative demonstrate that pseudo-BLUP should be extended to include collateral information sources beyond full-and half-sibs and their parents' EBV's.

Implications for selection schemes
The developed tool is a useful tool in assessing a short-term outcome of a selection scheme for a maternally affected trait. We considered simple examples on selecting for early growth performance. Using the parameters of Nordberg et al. (2005) for daily gain until two months' age (DG2) in sheep, let , , and the proportion of common environmental variance be 0.3, 0.1, -0.4 and 0.0 for DG2, respectively. Asymptotic response to selection were predicted with 25 ewes per ram and 3 lambs per (Finnsheep) ewe. In males and females, the proportion selected was 10 and 50%. The asymptotic response in DGE was 70-75% of the first generation response. With given genetic parameters, very high (1:3) economic weight was needed before any positive response was achieved in the maternal effect. In this case the asymptotic response in the direct effect was 13-29% lower than the one with equal economic weights. In reindeer research, our aim is to find an optimum data collection strategy for understanding the selection on average daily gain (ADG) from birth to autumn (cf Rönnegård et al. 2003). For the trait, Muuttoranta et al. (2014) found = 0.27 -0.52, = 0.22 -0.31 and = -0.74 --0.80. In reindeer there are in practice no full-sibs within a year class and the selection has to be based on calf's own record and those of its maternal half-sibs. Using = 0.40, = 0.25 and = -0.50, with equal weighting of the effects we need information on more (3-5) half-sibs, to have a positive response in MGE or high (1:3) economic weight. The reduction in the heritabilities and genetic correlation was almost negligible in both examples. Table 5. The asymptotic rate of response to selection (mean from 7 th to 10 th generation, in units of the base generation). There were 50 sires with f = 10 and n = 10. Selection intensity was 2.06 for both sexes. (A) The impact of changing the economic weights (w) for direct ( ) and maternal ( ) effect when = 0.3, and = 0.1, 0.3 or 0.5 and = -0.5. (B) The impact of common environmental variance ( ) when = = 0.3 and = 0.5 or -0.5. ℎ 2 = ℎ 2 = 0.3 ℎ 2 ℎ 2 = ℎ 2 = 0.3 ℎ 2 ℎ 2 = ℎ 2 = 0.3 ℎ 2 ℎ 2 = ℎ 2 = 0.3 ℎ 2 = ℎ 2 = 0.3 ℎ 2 = ℎ 2 = 0.3 ℎ 2 = ℎ 2 = 0.3 2 2 In reindeer, due to intermediate accuracy and delays over several generations in the evaluation of maternal genetic effects, an early selection for their improvement would benefit from the use of genomic information in parent diagnostics and ultimately in genomic selection. Lillehammer et al. (2020) used simulation in analyzing genomic selection for maternal traits in sheep with extant large scale recording scheme. In reindeer for these purposes, we need to develop a cheap SNP panel based possibly on the recent genome sequence information (Weldenegodguad et al. 2020). The marker information can be incorporated in the pseudo-BLUP as a correlated information (Dekkers 2007) which would provide a useful method to assess also the benefits of genomic selection and an optimum recording scheme (Mulder 2016).

Conclusions
The pseudo-BLUP for maternal traits provides a basis to assess the requirements for data collection and selection strategies. It is obvious that to handle negative genetic correlation between direct and maternal effects in approximating BLUP selection we need also more distant collateral relatives. The pseudo-BLUP is very flexible and can be extended to multi-trait situation and therefore for the consideration of several economically important traits and the utilization of genomic information.