Response to mass selection when the genotype by environment interaction is modelled as a linear reaction norm

A breeding goal accounting for the effects of genotype by environment interaction (G × E) has to define not only traits but also the environment in which those traits are to be improved. The aim of this study was to predict the selection response in the coefficients of a linear reaction norm, and response in average phenotypic value in any environment, when mass selection is applied to a trait where G × E is modelled as a linear reaction norm. The optimum environment in which to test the selection candidates for a given breeding objective was derived. Optimisation of the selection environment can be used as a means to either maximise genetic progress in a certain response environment, to keep the change in environmental sensitivity at a desired rate, or to reduce the proportion of animals performing below an acceptance level. The results showed that the optimum selection environment is not always equal to the environment in which the response is to be realised, but depends on the degree of G × E (determined by the ratio of variances in slope and level of a linear reaction norm), the correlation between level and slope, and the heritability of the trait.


INTRODUCTION
Genotype by environment interaction (G × E) is becoming increasingly important due to the globalisation of animal breeding. With G × E, the phenotypic expression of a trait in different environments, such as countries, climatic zones or production systems, is genetically not the same trait. In such cases, the breeding goal should define not only the traits but also the environment in which those traits are to be improved. The reaction norm model, where the phenotype is described as a continuous function of an environmental variable [14] is useful for studying G × E, especially when phenotypes change gradually over an environmental scale [8], e.g., production level or a climatic variable. In the range of environments normally encountered by a population of domestic animals, it is often reasonable to assume that reaction norms are linear functions of the environment, as has been found for milk production traits and fertility in dairy cattle [2,11].
In animal breeding, a substantial research effort has been devoted to data analysis using reaction norm models or statistically similar random regression models [2,7,11], whereas a limited effort has been made for the optimisation of breeding programmes for those situations. Kirkpatrick and Bataillon [10] have derived equations for the maximisation of selection response in the phenotypic value in a specified environment to mass selection on a trait affected by G × E. Their approach was to derive optimum index weights for observations recorded in different environments, modelling a covariance function without any assumptions of the shape of the reaction norm. Another approach is to model G × E as a linear reaction norm. The advantage of this model is that selection response can be predicted not only in the phenotypic expression in any environment, but also in the environmental sensitivity of the trait (robustness or responsiveness to changes in the environment, the slope of a linear reaction norm). Selection response in reaction norm coefficients has previously been described in terms of a selection gradient, expressing the covariance between the coefficients and fitness [8,9].
The objective of this study was to describe selection response of a trait affected by G × E in terms of the selection index theory. Equations were derived for the prediction of genetic change in reaction norm coefficients depending on the environment in which the animals were tested. The response in the average phenotypic value in any environment can be calculated knowing the genetic change in reaction norm coefficients. With mass selection and G × E, the variables available for the optimisation of selection response are selection intensity and the selection environment, which affect the accuracy of selection. The focus of this study was to derive prediction equations to find the optimum environment in which to test the selection candidates for a given breeding objective.

METHODS AND RESULTS
The selection index theory was combined with a reaction norm model and the Bulmer effect. Equations will be derived for the prediction of genetic change in reaction norm parameters at equilibrium genetic parameters.
From these equations, other equations will be derived in order to find the optimum selection environment for (i) maximising genetic progress in a certain environment, (ii) keeping the change in environmental sensitivity at a desired rate, or (iii) reducing the proportion of animals performing below an acceptable level. The equations were the main results of this study. The results are illustrated in connection to their derivation in the methods and results section. The implications of the results will be discussed in the discussion section.

The linear reaction norm model
The model was a linear reaction norm function for a single trait. In the following, the intercept and linear coefficient of a linear reaction norm will be referred to as level and slope. The phenotype of an individual j in an environment k was modelled as where y jk is the phenotypic value, b 0 is the population average level (intercept) in the average environment, b 1 is the population average slope, x k is the effect of environment k on the phenotype, a 0 j and a 1 j are the true breeding values for level and slope, respectively, and e 0 j and e 1 j are the environmental (residual) effects on level and slope, respectively. The term b 0 + b 1 x k represent the population average reaction norm. The intercept b 0 is positioned in the average environment so that E(x) = 0, and breeding and environmental values are expressed as deviations from the average reaction norm. As is common in animal breeding, covariances between residuals of different individuals and covariances between breeding values and residuals were assumed to be zero. If there are reasons to assume that this is not true, the model could be extended by the inclusion of an effect of a common environment, e.g., for individuals of the same litter. The slope of a linear reaction norm is a measure of sensitivity towards environmental change, which can be treated as a trait of the animal [8]. Genetic variation for the trait environmental sensitivity results in G × E and a genetic correlation <1 between phenotypic values of another trait measured in two different environments.
Note that both the genetic and the environmental effects are assumed to be linear functions of the environmental value. Consequently, also genetic and environmental variances change with the environment. The phenotypic variance in environment k, σ 2 y k , i.e., the variance of equation (1), was where x k is a row vector 1 x k of the environment k and G, E and P are the genetic, environmental and phenotypic (co)variance matrices, respectively, of the reaction norm parameters.
where σ 2 a 0 and σ 2 e 0 are the genetic and environmental variances of level, σ 2 a 1 and σ 2 e 1 are the genetic and environmental variances of slope, and σ a 0 a 1 and σ e 0 e 1 are the genetic and environmental covariances between level and slope, and P = G + E. The heritability in environment k becomes x k Gx k /x k Px k .

Population parameters
To illustrate the theoretical results that will be derived in the following sections, a Fortran 90 deterministic simulation programme was written. The infinitesimal model was assumed. Input values for the simulation were genetic and environmental parameters for the base population. The base population total phenotypic variance was set to 1.0 in the environment of the intercept of the reaction norm. Due to the variance of the slope, the genetic, environmental and, consequently, also the phenotypic total variance changed with the environment. The genetic and phenotypic total variances in three environments are shown in Table I. When there was no correlation between level and slope, the variance increased symmetrically with increasing distance from the environment of the intercept. With a non zero correlation between level and slope, the variance changed asymmetrically at the two sides of the environment of the intercept.
To illustrate varying degrees of G × E, the genetic correlation between the trait expressed in the environment of the intercept and an environment deviating 1 SD was set to 0.95, 0.80, and 0.60, corresponding to ratios of the genetic variances of slope and level of 0.11, 0.56, and 1.77, respectively, assuming the genetic and environmental correlations between level and slope are zero (Tab. I). Two other values, −0.4 and 0.4, were also studied for the correlation between level and slope. The correlation between level and slope had little Table I. The correlation between the level and slope (r (level, slope), genetic and environmental correlations assumed equal), variance ratio between slope and level, genetic correlation between the expression of a trait in the environment of the intercept and an environment deviating 1 SD from the environment of the intercept (r g (x = 0, x = 1)), and the base population total genetic variance (σ 2 G , assuming h 2 = 0.5 in the environment of the intercept) and total phenotypic variance (σ 2 P ) in the environment of the intercept (x = 0) and an environment deviating 1 SD from the environment of the intercept (x = −1 and x = 1) assuming h 2 was constant over the range of environments. influence on the genetic correlation between the trait expressed in the environment of the intercept and an environment deviating 1 SD when the variance in the slope was small (a small variance ratio, little G × E), but a larger influence with a larger variance in the slope (larger variance ratio, more G × E) (Tab. I). The genetic and environmental correlations between level and slope were assumed to be equal. Two values of base population heritability in the environment of the intercept were studied: 0.2 and 0.5. For the illustration of selection response in level and slope and in Table I, the genetic and environmental variance ratios of slope and level were set to be equal. Then the genetic and environmental variances were affected proportionally by the environmental effect and heritability was constant over the environmental range. For the illustration of selection response of the total phenotypic value in a specified environment and of the optimum selection environment, the environmental variance of slope was increased so that heritability decreased with an increasing distance from the intercept, being 5% lower at 1 SD from the intercept.
The Bulmer effect was accounted for by an iterative reduction of the genetic parameters until equilibrium was established (App. A). The equilibrium values were used for calculations of selection response, accuracy of selection, and optimum selection environment. Mass selection with 10% of the males and 10% of the females selected was assumed. The selection response is expressed in phenotypic SD-units of the trait expressed in the environment of the intercept.

Genetic change in reaction norm parameters
The genetic change in the reaction norm parameters, level and slope, when selecting on the phenotypic value y k is a function of the selection environment, k. From the regression of the breeding values for level and slope on the phenotypic selection differential in environment k, it follows that where ∆a 0 and ∆a 1 are the selection responses in level and slope, respectively, i is the selection intensity, and σ y k = x k Px k is the phenotypic standard deviation in environment k.
The accuracy of selection is the correlation between the selection criterion (the phenotypic value) and the reaction norm parameter and for level and slope, respectively, and where r a 0 ,a 1 is the genetic correlation between level and slope. The accuracies can be substituted into equation (3) to give ∆a 0 = ir y k ,a 0 σ a 0 and ∆a 1 = ir y k ,a 1 σ a 1 , which agrees with the classical ∆G = ir IA σ A where r IA is the correlation between the index and the true breeding value and σ A is the genetic standard deviation [5].
The selection responses in level and slope over a range of environments are shown in Figure 1 for a heritability that is constant over the environmental scale. As expected, the selection environment had a larger effect on the selection response with a higher variance in slope, i.e., a the higher degree of G × E. The correlation between level and slope affected the shape and location on the environmental axis of the response curve. The maximum selection response in level was achieved with selection in the average environment. When the correlation between level and slope was positive, the response in level was higher with a selection environment that was better than the average environment than in a selection environment that was worse, and vice versa with a negative correlation. Selection response in slope was always the highest with a selection environment that was better than average and with a positive correlation between level and slope. With the lower level of heritability (0.2) the response in both level and slope was smaller and less affected by the environmental value, but the shape of the response curve was the same as with the higher level of heritability (0.5) (not shown). The intermediate level of G × E (variance ratio 0.56, not shown) yielded results that were intermediate between the high and low level of G × E shown in Figure 1.

Genetic change in other environments
The genetic change in a defined environment, l, is a function of both the environment of selection, k, and the environment, l, where the results of selection are expressed. When multiplying equation (3) by x l we get where ∆G l is the genetic change in environment l, x l is a column vector 1 x l , and x k Gx l is the covariance between the selection criterion y k and the genetic merit in environment l. The accuracy of selection, r kl , is the correlation between the selection criterion y k and the genetic merit in environment l, a 0 j + a 1 j x l .
where σ A l = x l Gx l is the standard deviation of the genetic merit in environment l. Note that combining equations (6) and (7) gives ∆G l = ir kl σ A l , as expected.
The genetic correlation, r g 1,2 , between the genetic merit in two environments 1 and 2, is where x 1 Gx 2 is the covariance between the genetic merit in environment 1 and 2, and x 1 Gx 1 and x 2 Gx 2 are the standard deviations of the true breeding values in environments 1 and 2, respectively. The selection response in environment l as a function of the selection environment k is illustrated in Figure 2 for a heritability of 0.5 in the environment of the intercept and 5% lower 1 SD from the intercept. Figure 2 shows that maximum gain was achieved when the selection environment was close to the response environment, but not necessarily equal. The correlation between level and slope affected the shape of the response curve. Heritability (0.2 or 0.5 in the environment of the intercept) affected the magnitude of response (not shown). As for the response in level and slope, the selection environment had a larger effect on the selection response when there was a higher variance in slope / and a higher G × E (environmental effect with variance ratio 1.77 > 0.56 > 0.11, variance ratio 0.56 not shown). This sensitivity of the response to the selection environment was asymmetric when level and slope were correlated, i.e., the cost of a sub-optimal selection environment depended on which side of the optimum environment was the actual selection environment. With a constant heritability the maximum genetic gain was achieved when the selection environment was equal to the response environment (not shown).

Maximising genetic progress
A breeding goal that defines a goal trait and the environment in which the trait is to be improved can be expressed as the expected phenotypic value, y jl , in environment l where the animals are expected to perform in the future, E(y jl ) = a 0 j + a 1 j x l . The selection criterion is the phenotypic value, y jk , in the selection environment k. The optimum selection environment, k opt , that gives the highest ∆G l is found by solving ∂∆G l /∂x k = 0. The derivation is shown in Appendix B.
For the general case, When there is neither a genetic nor an environmental correlation between level and slope the expression simplifies to x k opt = x l , when h 2 0 is defined as σ 2 a 0 /σ 2 p 0 , and h 2 1 as σ 2 a 1 /σ 2 p 1 . When the genetic and environmental variance ratios are equal so that heritability is constant over the environmental range x k opt simplifies to x k opt = x l . Note that the latter result agrees with the classical expression for correlated response, CR l = ih k r g k,l σ A l [6]. At a constant heritability, the response in trait expression in environment l, CR l , to selection on trait expression in environment k is maximised when the genetic correlation between trait expression in the two environments, r g k,l is unity, i.e., the environments are equal.
Equation (9) shows that the optimum selection environment is a function of the equilibrium genetic parameters. The equilibrium genetic parameters depend on the efficiency of selection which, in turn, depends on the selection environment. Therefore the optimum selection environment depends on the environment in which the equilibrium genetic parameters have been obtained. To find the equilibrium optimum selection environment, an optimum selection environment was first calculated using equation (9) and base population parameters. Then the equilibrium genetic parameters in this environment were calculated and used to calculate a new optimum environment. This procedure was iterated until the optimum environment did not change between rounds. In our simulation, equilibrium was reached after 1−5 iterations with a convergence criterion of less than 0.1% change in an optimum environment between rounds.
The equilibrium optimum selection environment for different response environments is shown in Figure 3a. The optimum selection environment was not equal to the response environment and their relationship was not precisely linear. With more G×E (higher variance ratio) and/or a flatter heritability function (2% difference between heritability in average and deviating environment, not shown) the optimum selection environment was closer to the response environment, than with the lower variance ratio and/or a more curved heritability function.
To investigate the need for an iterative solution of x k opt , x k opt was calculated also using base generation parameters, and the results were compared to the solution using equilibrium parameters. Figure 3b shows that the optimum selection environment calculated from base population parameters was close to the equilibrium optimum selection environment. With a higher variance in slope, the difference between the optimum selection environment calculated using the base population or equilibrium parameters was even smaller than that shown in Figure 3b. Thus the calculation of x k opt directly from base generation parameters is probably sufficient for practical applications.

Keeping environmental sensitivity constant
Selection for high phenotypic value in combination with a continuous improvement of the environment and in the presence of G × E has been shown to increase the environmental sensitivity of a population [12]. In modern animal production systems, genetic improvement and improvements of feeding and management occur simultaneously. A restriction of the genetic change in environmental sensitivity can be included in the breeding objective to avoid that the animals become increasingly sensitive to, e.g., feed quality or disease.
Since the genetic change in slope is a function of the selection environment, k, it is possible to find a selection environment that gives no change in slope. This environment can be found numerically from a plot of ∆a 1 over x k (Fig. 1b). It can also be found analytically as the environment where there is no covariance between the selection criterion and the true breeding value for the slope. From equation (5) we get σ y jk a 1 = σ a 0 a 1 + x k σ 2 a 1 . Equating this to zero gives x k ∆a 1 =0 = −σ a 0 a 1 /σ 2 a 1 or, rearranged, x k ∆a 1 =0 = −r a 0, a 1 σ a 0 /σ a 1 . In words, the selection environment that gives no change in slope is "minus the correlation between level and slope, scaled by the ratio of the genetic standard deviations of level and slope". Note that, if there is a positive correlation between level and slope, x k ∆a 1 =0 must be negative, i.e., worse than the average environment.

Desired gains in level and slope
Keeping the slope constant, as in 2.3.4., is a special case of desired gains for level and/or slope. For the general case, a desired ratio, δ, between the change in slope and level can be specified. From equation (3), it follows that Solving for x k δ , the selection environment where this  ratio is achieved, gives x k δ = (σ a 0 a 1 − δσ 2 a 0 )/(δσ a 0 a 1 − σ 2 a 1 ). If no change in slope is desired, then δ = 0 and x k δ = −σ a 0 a 1 /σ 2 a 1 , as shown in the previous section.

Reducing the proportion of animals performing below an acceptable level
The main interest might not be in all cases to increase the level or to reduce the slope. For traits related to animal welfare, it may be more relevant to avoid poor animal performance. Besides the suffering of the animal and potential costs and worries for the farmer, poor animal performance may reduce consumer acceptance of the production system. Therefore, an objective may be to reduce the proportion, p, of animals that perform below an acceptable phenotypic level, say τ y (Fig. 4).
Consider the population average reaction norm, E(y k ) = b 0 + b 1 x k . The animals in an environment x k < τ x = (τ y − b 0 )/b 1 are expected to perform below the acceptance level τ y , where τ x is the environmental truncation point, below which y k|τ x falls below the phenotypic threshold τ y = b 0 + b 1 τ x . The breeding objective, H, is to reduce the proportion p, thus H = −p, which is a non-linear function of level and slope. The selection criterion is the phenotypic value in the selection environment k. The problem is to find the selection environment that yields the largest reduction of p. The first step is to find the optimum weights, v 0 and v 1 , for level and slope. The breeding objective can be approximated by a linear function of level and slope, H ≈ v 0 a 0 + v 1 a 1 . Values for v 0 and v 1 are obtained by taking the partial derivatives of −p with respect to level and slope, and calculating the values of the derivatives at the current population mean (b 0 and b 1 ). The relationship between reaction norm parameters and p arises solely via τ x , so that partial derivatives can be obtained via τ The partial derivatives are ∂(−p)/∂τ x = −ϕ(τ x ), ∂τ x /∂b 0 = −1/b 1 , and ∂τ x /∂b 1 = (b 0 − τ y )/b 2 1 , where ϕ is the density function of the distribution of the environmental values. The resulting weights are v 0 = ϕ(τ x )/b 1 Since the scaling of the weights does not affect the selection response, v 0 and v 1 may be divided by ϕ(τ x )/b 1 to give v 0 = 1 and v 1 = (τ y −b 0 )/b 1 = τ x . Thus, the linearised breeding objective is H = a 0 +a 1 τ x . Note that this breeding goal corresponds to genetic gain in environment τ x , which shows that maximising the reduction of p is the same as maximising the selection response in environment τ x . Consequently, the environment of selection that maximises −∆p can be found by substituting x l by τ x in equation (9).
The result that minimising p is the same as maximising the genetic change in environment τ x can be derived also in a more intuitive manner as follows: The intercept of a reaction norm can be defined in any environment without changing the reaction norm itself (see the discussion section). If the intercept is defined at the environmental truncation point τ x , instead of in the average environment, the rotation of the reaction norm (change in slope) around this intercept does not affect τ x . For small changes in level and slope, the shift of the truncation point τ x , and hence, the change in p is determined entirely by the change in level/phenotypic value in the environment τ x . This argument extends directly to higher-order reaction norms, meaning that also for quadratic and cubic reaction norms, maximising −∆p is identical to maximising ∆G in the environment τ x . Equation (9) will need extension to higher-order reaction norms, though.

DISCUSSION
The described theory is a tool for understanding the effects of G × E on selection response as a function of the environment of selection and the environment where selection response is to be realised. In the case of mass selection, the derived equations can be used to solve the problem stated by Falconer [4] whether to select under the conditions in which the breed is required to live, or if it would be better to select "under some other conditions, for example under conditions more favourable for the expression of the desired character". The old thesis that "the character required is best selected for under environmental conditions which favour its fullest expression" [6] may suggest that one should select in the environment, which maximises heritability. Equation (9), however, shows that the optimum selection environment depends on, besides the environment where selection response is to be realised, three factors: the degree of G × E (determined by the ratio of variances in slope and level), the correlation between level and slope, and the heritability of the trait. The results show that the optimum selection environment is neither always equal to the environment in which selection response is to be realised, nor to the environment where heritability is the highest (Fig. 2).
The numerical illustrations show the relations between genetic and environmental parameters and selection response. The ratio between the variances in slope and level determines the magnitude of the environmental effect on selection response. With a high variation in slope, i.e., a high degree of G × E, the environmental value has a large effect on selection response. For prediction of selection response for a few generations of selection, which is often of interest in animal breeding, the prediction equations derived in this paper are relevant even when the conditions of the infinitesimal model do not hold exactly.
The equations presented here enable the breeder to quantify the value of data recording in a particular environment, and also to quantify the losses due to data recording in sub-optimal environments. The small difference between the optimum selection environment calculated using a base population or equilibrium parameters suggests that the optimum environment can, for convenience, be calculated from base population parameters without losing much in genetic gain.
It is important to realise that with a reaction norm the correlation between level and slope depends on the location of the intercept on the environmental scale [13]. Consider the phenotypic reaction norm y jk = b 0 +b 1 x k + p 0 j + p 1 j x k , where p 0 j and p 1 j are the phenotypic values of level and slope of individual j. A non-zero correlation between p 0 j and p 1 j can be removed by redefining the intercept in such a way that the phenotype y jk at the new intercept has no covariance with p 1 j . The environment where this occurs is found by solving σ y jk ,p 1 j = σ p 0 j +p 1 j x k ,p 1 j = 0, which gives x k = −σ p 0 p 1 /σ 2 p 1 . If the original intercept is at x k = 0, then the environment where level and slope are not correlated is at x k = −σ p 0 p 1 /σ 2 p 1 . Hence, the new value of the intercept is b * 0 = b 0 − b 1 σ p 0 p 1 /σ 2 p 1 , and the slope is the same as for the original reaction norm. Thus, to remove a correlation between level and slope, a positive correlation requires shifting the intercept an amount σ p 0 p 1 /σ 2 p 1 downwards on the environmental scale, whereas a negative correlation requires shifting the intercept the same amount upwards on the environmental scale. A zero correlation can be achieved simultaneously on the genetic and environmental levels only when r a 0 a 1 = r e 0 e 1 and heritability is constant over the environmental scale.
Due to the dependency of the correlation between level and slope on the location of the intercept on the environmental scale, estimates of the correlation between level and slope should be interpreted with caution. To facilitate the comparison of results from different reaction norm studies, authors should report where the intercept was defined. We agree with the recommendation of van Tienderen and Koelewijn [13] to define the intercept in the average environment.
For species where the breeding population comprises almost the entire population, such as dairy cattle, the selection environment cannot simply be chosen. With progeny testing, the accuracy is generally rather high, so that the optimum environment for data recording will depend less on heritability and more on the genetic correlation with the breeding goal. Equation (9) does not apply to progeny testing. Instead one can argue that the optimum environment to record progeny information will be close to the environment where the population has to perform in the future. Then, when progeny testing of sires is combined with a nucleus herd for selection of cows, a near optimum solution would be to have the nucleus animals in an environment as indicated by equation (9), and the progeny in an environment similar to the expected future environment of production animals.
A problem arises if the environment that would give the maximum selection response is not an ethically acceptable environment for testing animals. Testing animals in the optimum environment for improving disease resistance, or to reduce the proportion of animals performing very poor, might be comparable to doing a challenge test. A major benefit in animal health and farmer economy in coming generations could justify such a test, though for practical as well as ethical reasons it should be done only with a small part of the population.
In order to be more applicable to practical animal breeding, the theory presented here needs to be extended to include multiple trait evaluations, more complicated selection schemes, such as selection on BLUP procedure breeding values, and non-linear reaction norm functions. Here we have assumed that animals are in a single environment, whereas a distribution of animals across environments is more realistic. Such a model would be useful for searching for