## Abstract

The differences of a phenotypic trait produced by a genotype in response to changes in the environment are referred to as phenotypic plasticity. Despite its importance in the maintenance of genetic diversity via genotype-by-environment interactions, little is known about the detailed genetic architecture of this phenomenon, thus limiting our ability to predict the pattern and process of microevolutionary responses to changing environments. In this article, we develop a statistical model for mapping quantitative trait loci (QTL) that control the phenotypic plasticity of a complex trait through differentiated expressions of pleiotropic QTL in different environments. In particular, our model focuses on count traits that represent an important aspect of biological systems, controlled by a network of multiple genes and environmental factors. The model was derived within a multivariate mixture model framework in which QTL genotype-specific mixture components are modeled by a multivariate Poisson distribution for a count trait expressed in multiple clonal replicates. A two-stage hierarchic EM algorithm is implemented to obtain the maximum-likelihood estimates of the Poisson parameters that specify environment-specific genetic effects of a QTL and residual errors. By approximating the number of sylleptic branches on the main stems of poplar hybrids by a Poisson distribution, the new model was applied to map QTL that contribute to the phenotypic plasticity of a count trait. The statistical behavior of the model and its utilization were investigated through simulation studies that mimic the poplar example used. This model will provide insights into how genomes and environments interact to determine the phenotypes of complex count traits.

ONE of the most important challenges facing modern biology is to understand the genetic mechanisms underlying the adaptation of biological traits to environmental factors and use this knowledge to predict the response of biological structure, organization, and function to changing environments (Franks *et al*. 2007; Huntley 2007). When grown in different environments, an organism may show a range of phenotypes. Such a capacity of the organism to alter its phenotypes in response to changing environment, defined as phenotypic plasticity, has been recognized by many early biologists (Waddington 1942; Schmalhausen 1949). Some earlier theoretical models, including homeostasis, were useful for explaining why some individuals are insensitive to changes in environment (Waddington 1942). With a rapid and dramatic change in global climate stimulated by human activities (Grether 2005), the study of phenotypic plasticity has received a resurgence of interest and a renaissance in its fundamental role in shaping evolutionary adaptation and consequences (Scheiner 1993; Via *et al*. 1995; Schlichting and Smith 2002; West-Eberhard 2003, 2005; Wu *et al*. 2004; De Jong 2005). Because different genotypes display a wide range of variation in the level of their plasticity response (Ungerer *et al*. 2003), there must be a genetic basis for phenotypic plasticity to environmental change. For this reason, the identification of specific genes that contribute to plastic responses has now become an important research area for understanding the genetic and developmental machineries of organismic adaptation and evolution (Gibert *et al*. 2007).

Genetic mapping of complex traits with molecular markers has been proven powerful for the genomewide characterization of quantitative trait loci (QTL) that regulate phenotypic plasticity (Wu 1998; Leips and Mackay 2000; Kliebenstein *et al*. 2002; Ungerer *et al*. 2003; Gutteling *et al*. 2007). The statistical test that measures differences in genetic effects of a QTL across different environments provides meaningful procedures for investigating the genetic mechanisms hypothesized to explain the genetic basis of phenotypic plasticity. The overdominance hypothesis proposes that a heterozygote at relevant genes shows higher stability, or lower plasticity, than a homozygote and the degree of stability is proportional to heterozygosity for these loci (homeostasis; Gillespie and Turelli 1989). The pleiotropic hypothesis states that differential expression of the same loci across environments causes phenotypic plasticity (allelic sensitivity; Via and Lande 1985). The epistatic hypothesis suggests that specific plasticity genes exist that interact epistatically with the loci for the mean value of the trait to regulate environmental sensitivity (gene regulation; Scheiner and Lyman 1989). These hypotheses have been tested with results from QTL mapping in different species from Arabidopsis (Kliebenstein *et al*. 2002; Ungerer *et al*. 2003) to Populus (Wu 1998; Rae *et al*. 2008), Drosophila (Leips and Mackay 2000; Geiger-Thornsberry and Mackay 2002; Anholt and Mackay 2004), and Caenorhabditis (Gutteling *et al*. 2007). All these studies were based on the phenotypic plasticity of continuously varying traits. The genetic control of phenotypic plasticity for count traits—another group of important traits to agriculture, biology, and biomedicine—is still poorly understood.

Statistical models for genetic mapping of continuous traits that are normally distributed have been well developed in the past two decades (Lander and Botstein 1989; Jansen and Stam 1994; Zeng 1994; Kao *et al*. 1999; Wu *et al*. 2007; Rae *et al*. 2008). The idea of mapping continuous traits has been extended to map binary or ordinal traits that vary in a discontinuous manner on the basis of a threshold model by assuming a continuously distributed liability underlying these traits (Visscher *et al*. 1996; Xu and Atchley 1996; Yi and Xu 1999; Li *et al*. 2006). Different from continuous traits (taking an infinite number of trait values) and ordinal traits (taking a fixed number of discrete trait values), there is also a group of traits measured in counts that are discrete yet may take an infinite number of values. Count data, such as cell numbers, branch numbers, or bristle numbers, play a unique role in determining the phenotypic plasticity of a biological system (Norga *et al*. 2003).

More recently, Cui *et al*. (2006) generalized a parametric mapping strategy, as proposed by Rebaï (1997), Shepel *et al*. (1998), and Sen and Churchill (2001), to map and test QTL controlling count traits. Different from traditional treatments based on normality or threshold assumptions, Cui *et al*. (2006) and others (Rebaï 1997; Shepel *et al*. 1998; Sen and Churchill 2001) modeled count data by incorporating an intrinsic Poisson distribution, thus leading to more reasonable biological interpretations about the genetic effects of a “count” QTL. Cui *et al*.'s model can be used as a standard procedure for mapping QTL contributing to the genetic control of complex count traits.

In this article, we integrate the Poisson distribution into a general mapping framework for the identification of QTL that control the phenotypic plasticity of a count trait through their environment-dependent expression. A particular study of phenotypic plasticity should make use of the advantages of a randomized complete block design with two or more different treatment levels and multiple replicates per each treatment level. Traditional interval mapping generally takes the means of a phenotypic trait over different replicates under each treatment level and then compares the difference in the genetic effect of each detected QTL across different treatment levels (Hayes *et al*. 1993). However, by taking averages over replicates, the averaged count traits may no longer be integer valued. Additionally, the covariance structure among the replicates cannot be utilized when averages are taken. To overcome these drawbacks, count observations in individual replicates are incorporated into the mapping model by invoking a multivariate Poisson distribution with dimension equal to the number of replicates. The estimation of the parameters that describe the multivariate Poisson distribution is obtained by the EM algorithm. The implementation of this algorithm into a mixture-based mapping model generates a two-stage hierarchical EM algorithm in which the Poisson parameters that define the environment-dependent genetic effects of a QTL can be estimated. The new model was used to map the QTL that affect sylleptic branch counts and their plasticity across two different fertilization regimes in a Populus hybrid population. Computer simulations are further used to study the performance of the method for mapping environment-sensitive QTL for complex count traits.

## MODEL

#### Model structure and estimation:

Consider a backcross population with *n* progeny in which there are two different genotypes at each locus. We assume that a genetic linkage map is constructed with polymorphic markers for this backcross. In many plants, such as poplar trees, clonal propagation is possible, thus allowing the same progeny to be genotypically replicated. Suppose the backcross considered is planted in a randomized complete block design with two treatment levels (*e.g*., low and high fertilization) and *R* clonal replicates within each treatment level. Each of the plants studied is measured for a count trait of interest, *e.g*., the number of branches on a main stem. This experimental design allows the characterization of genetic factors that control the response of each backcross progeny to different treatment levels.

Suppose there is a putative QTL segregating with two different genotypes *Qq* (coded by **1**) and *qq* (coded by **2**) in the backcross that affects the phenotypic plasticity of the trait across two treatment levels. This QTL is located somewhere in the genome, which can be detected by the linkage map. Assume the QTL resides between a pair of flanking markers **M**_{1} and **M**_{2} each with two genotypes coded by 1 and 0. For each backcross progeny, it may carry one (and only one) QTL genotype, **1** or **2**. The probability with which a particular progeny (*i*) carries QTL genotype **1** or **2** depends on the marker genotypes of this progeny at the two flanking markers (**M**_{1} and **M**_{2}) that bracket the QTL. Under the assumption of independent crossovers, the probability of a QTL genotype given a marker genotype can be derived in terms of the recombination fractions between **M**_{1} and QTL, between QTL and **M**_{2}, and between the two markers. Given that each progeny has a known marker genotype, 11, 10, 01, or 00, the conditional probability of QTL genotype *Qq* for a given progeny *i* given its marker genotype is denoted as ω_{1}_{|i} and the conditional probability of QTL genotype *qq* is ω_{2}_{|i} = 1 − ω_{1}_{|i} (Wu *et al*. 2007).

The trait values of backcross progeny *i* in *R* different replicates under treatment level *k* (*k* = 1, 2), arrayed in **X**_{ik} = (*X _{ik}*

_{1}, … ,

*X*), are distributed as a mixture function with two different groups of QTL genotypes;

_{ikR}*i.e*.,(1)where

**Θ**

_{k}= (

**Θ**

_{1}_{|k},

**Θ**

_{2}_{|k}) contains parameters specific to QTL genotype

*j*for treatment

*k*, and

*P*(

_{j}**X**

_{ik}|

**Θ**

_{j|k}) is a probability density function for the count trait, which can be described by a multivariate Poisson distribution, expressed as(2)where each

*X*follows a Poisson distribution with mean parameter θ

_{ikr}_{j|kr}that is the genotypic mean of the count trait for QTL genotype

*j*in replicate

*r*at treatment level

*k*and covariance parameter θ

_{j|k0}that is the QTL genotype- and treatment level-specific covariance of the count trait between all the pairs of replicates. If θ

_{j|k0}= 0, then the variables are independent and the multivariate Poisson distribution reduces to the product of independent Poisson distributions.

Assuming that the trait values from different levels of treatment are independent, we construct a likelihood of the unknown parameters given the trait values and marker information (**M**) in terms of the mixture model (1) by combining the two treatment levels; that is(3)

The maximum-likelihood estimates (MLEs) of the parameters can be obtained by maximizing the likelihood (3). In this clonal design, it is reasonable to assume that genotypic means of the count trait are equal among different replicates at each treatment level; that is,(4)Equation 4 suggests that, unlike the general multivariate Poisson model in which different variables are assumed to have different means, our QTL mapping model assumes the same mean for all the variables. Thus, the set of parameters being estimated in the likelihood (3) is **Θ**_{j|k} = (θ_{j|k0}, θ_{j|k}) (*j* = **1**, **2**; *k* = 1, 2). The EM algorithm can be implemented to estimate the maximum-likelihood estimates (MLEs) of these parameters as follows.

##### E step:

Given the data and the current values of the estimates after the *t*th iteration and , we calculate the conditional expectations of the complete data, which include the pseudovalues of progeny *i* within treatment level *k*,and the posterior probabilities, with which progeny *i* at treatment level *k* carries QTL genotype *j*,

##### M step:

The estimates of parameters are updated by usingwhere . Both the E and the M steps are iterated until the estimates converge to stable values according to some prespecified convergence criterion. As usual, the QTL position is estimated via a fixed approach with which the existence of a QTL is scanned at every 2 cM within a marker interval. Significant peaks of the log-likelihood ratio value profile across the genome are thought to be the estimated location of a QTL (see below).

#### Hypothesis testing:

After the genetic parameters are obtained, we need to test whether there is a QTL that affects the count trait at the two levels of treatments. The existence of a QTL can be tested by formulating the hypotheses(5)where the null hypothesis H_{0} states that the data can be fit with only one mean for each treatment level, whereas in the alternative hypothesis H_{1} there are two distinct means showing that there is a segregating QTL for the trait. The test statistic is the log-likelihood ratio (LR) of the full (H_{1}) over reduced model (H_{0}), expressed aswhere the tildes and the circumflexes denote the MLEs of the unknown parameters under H_{0} and H_{1}, respectively. Note that the estimation of depends on both phenotypic values and marker data, whereas the estimation of depends only on phenotypic values. The critical threshold for the declaration of a QTL can be determined from permutation tests (Churchill and Doerge 1994).

If a significant QTL is found, then we can test whether this QTL has a different effect on trait values at two different treatment levels. This test is formulated as(6)and(7)If the null hypotheses above are both rejected, this implies that the significant QTL detected affects pleiotropically trait phenotypes at two different treatment levels. Otherwise, this QTL is operational only at one treatment level, which causes genotype-by-environment interactions for the trait.

In practice, although the QTL is pleiotropic, its effect on the trait may depend on the level of treatment. This means that the same QTL has different effects in different environments. Differential expression of a QTL across two different environments can be tested using(8)The rejection of the above null hypothesis means that a significant genotype-by-environment interaction exists due to allelic sensitivity to a varying environment. Significant genotype-by-environment interactions are one of the genetic causes of phenotypic plasticity. The log-likelihood-ratio test statistics for hypotheses (6)–(8) can be determined from simulation studies.

## RESULTS

#### A worked example:

An interspecific hybrid family (designated as 52-124) was generated by crossing *Populus deltoides* (genotype ILL101) to *P. trichocarpa* (genotype Gf93-968). One F_{1} female clone (52-225) was then crossed with a different *P. deltoides* individual (D124), which resulted in a pseudobackcross progeny of 397 individuals. Family 52-124 was genotyped with 171 microsatellite (SSR) markers. The linkage map was constructed using MapMaker version 3.0 (Lander *et al*. 1987) following the two-way pseudotestcross mapping strategy (Grattapaglia and Sederoff 1994), generating a separate map for each of the parents. Markers were assigned to linkage groups at a minimum threshold of LOD 3 and a recombination fraction of 0.35. The resulting linkage was mapped onto 19 linkage groups, equivalent to the Populus chromosome number, with all markers displaying internally consistent linkage patterns.

Greenwood cuttings (clonal replicates) from each individual of family 52-124 were produced. All plants were initially grown for 3 weeks under 5 mm of nitrogen (N) fertilizer, in a modified Hockings solution (Cooke *et al*. 2003). Next, plants were subjected to two different N treatment levels (0 and 25 mm) in a randomized complete block design, with three replicates per level. The number of sylleptic branches derived from the main stem (*i.e*., those derived from lateral buds without an intervening stage of rest; Wu and Hinckley 2001) were counted for each tree 30 days later.

We use the pseudotest backcross mapping strategy proposed by Grattapaglia and Sederoff (1994) to map branch numbers. This strategy constructs two different sets of linkage maps with molecular markers that are segregating in one parent but not in the other. Although it does not take into account all types of markers segregating in the family, the pseudotest backcross strategy has been useful for mapping QTL controlling quantitative traits in forest trees (Grattapaglia *et al*. 1996). The number of branches on the main stems was observed to approximately follow a Poisson distribution under each of the fertilization levels (Figure 1). For a possible use of traditional normality-based approaches, the data were subjected to log transformation, but the transformed data still tended to obey a Poisson distribution rather than a normal distribution. Thus, existing QTL mapping models based on normality assumption may not work in this particular case. There are many more branches when trees are planted under a high than a low fertilization level, but the degree of increase for branch number with higher fertilization level seems to vary among different progeny (Figure 2), suggesting that genotype-by-environment interactions are important for branch number. The new method was used to scan for the existence and distribution of QTL that control branch numbers throughout father- and mother-based linkage maps using a joint likelihood (3). A significant QTL was detected between markers G2431 and P2855 on the mother-based linkage map since the peak of the log-likelihood ratio profile (1493.5) is beyond the genomewide critical threshold (1391.1) determined from 1000 permutation tests (Figure 3).

This significant QTL detected was further tested for its effects on branch number separately under each fertilization level according to hypotheses (6) and (7). It displays a significant effect under high fertilization (*P* = 0.02), but it is marginally significant under low fertilization (*P* = 0.05). The additive genetic effect of this QTL estimated from the pseudotest backcross is 5.968 under high fertilization, compared with 2.777 under a low fertilization level (Table 1). A hypothesis test with Equation 8 suggests that these two additive effects are significantly different (*P* = 0.03), suggesting the existence of a QTL-by-environment interaction (Figure 4). This QTL explains 49.5 and 70.2% of the total phenotypic variance for branch number under low and high fertilization levels, respectively.

#### Monte Carlo simulation:

To examine the statistical properties of the new model, we performed simulation studies under different scenarios. Consider a backcross for which a linkage group of length 200 cM was simulated with 11 equally spaced markers. Two different sample sizes (*n* = 100 and 300) were assumed to simulate the backcross. The backcross progeny was planted in a complete randomized block design with two different treatment levels and three clonal replicates per treatment level. A QTL that determines a count trait was hypothesized at 30 cM from the first marker of the linkage group. The phenotypic values of each backcross at each treatment level were simulated by assuming a trivariate Poisson distribution with a mean depending on the genotype of the assumed QTL. The genetic effect of a QTL is expressed as the difference between two QTL genotypes, which is described by θ_{1}_{|k} − θ_{2}_{|k} at treatment level *k*. Three different simulation scenarios were considered: (1) The genetic effect of a QTL is large and there are covariances among the three replicates, (2) the genetic effect of a QTL is large and there is no covariance among the three replicates, and (3) the genetic effect of a QTL is small and there are covariances among the three replicates. Simulation under each scenario was repeated 1000 times to estimate the means and mean square errors of the parameter estimates.

Table 2 tabulates the results of QTL position and effect estimation by the new model. The model can provide reasonable estimates of all the parameters with a modest sample size (100) under the three simulation scenarios considered. The estimation accuracy and precision can increase dramatically when the sample size increases to 300. The Poisson parameters can be better estimated when there is no covariance among replicates (scenario 2) than when such a covariance exists (scenario 1), but the reverse is true for the estimation of the QTL position. As expected, when two QTL genotypes diverge more largely, the QTL position can be better estimated (scenario 1 *vs*. 3). The power to detect significant QTL-by-environment interactions under different simulation scenarios was also calculated. Such power is the largest under scenario 1, followed by scenarios 3 and 2. In all the cases, increasing sample sizes can improve the power for the detection of QTL-by-environment interactions.

We conducted an additional simulation study in which multiple QTL are involved in the genetic control of a count trait. We simulated a genome composed of four different chromosomes (A–D), each 200 cM long, covered by 11 equally spaced markers. Three QTL are located at 45, 95, and 155 cM from the first marker on the first three chromosomes, respectively. The phenotypic data were simulated as trivariate Poisson random variables based on different QTL genotypes and each given treatment level. Due to the assumption in Equation 4, the simulation for a given treatment level *k* requires only two parameters, the covariance and mean parameter , for a given three-QTL genotype *j*_{1}*j*_{2}*j*_{3} (*j*_{1} = 1, 2 for the first QTL on chromosome A, *j*_{2} = 1, 2 for the second QTL on chromosome B, and *j*_{3} = 1, 2 for the third QTL on chromosome C). Assuming no epistasis, any three-QTL genotypic mean can be expressed as the sum of genotypic means at individual QTL; that is,where *a*_{1}, *a*_{2}, and *a*_{3} are the additive genetic effects of the simulated QTL, respectively. Three QTL were assumed to display different effects at two different levels of treatment. Figure 5 illustrates the profile of the LR values calculated from an arbitrary simulation run to test the distribution of QTL across all four simulated chromosomes. As shown, all three QTL can be detected from our model. Table 3 lists the MLEs of the positions and genetic effects of different QTL and the square roots of mean square errors of the MLEs when different sample sizes are assumed. Our model can estimate the positions of three different QTL and their genetic effects with reasonable accuracy and precision. As expected, QTL of larger effects can be better estimated than those of smaller effects. Also, the estimation precision can be increased with increased sample sizes.

## DISCUSSION

Genetic mapping of QTL has been instrumental for the characterization of major loci that are responsible for a variety of quantitative traits (Mackay 2001; Anholt and Mackay 2004; Li *et al*. 2006; Paterson 2006). More recently, QTL mapping techniques have been considerably extended to approach many biologically meaningful issues, among which genetic mapping of the phenotypic plasticity of a count trait has not received adequate attention as compared to its fundamental importance in statistics and biology. While most statistical models for QTL mapping are based on the normality assumption for continuous traits (Lander and Botstein 1989; Zeng 1994; Jiang and Zeng 1995; Lynch and Walsh 1998; Wu *et al*. 2007) or the log-normality assumption for binary or ordinal traits (Li *et al*. 2006), there is also a group of traits in nature and of importance to agriculture and biomedicine that are expressed in counts, such as the numbers of leaves, branches, and seeds produced by a plant or the number of cancer cells. These count traits are typically better described by a Poisson distribution than by a normal distribution, and thus many existing models may not be appropriate for mapping this type of trait. To fully consider the statistical nature of count traits, some authors attempted to incorporate the Poisson distribution into a standard mapping model (Rebaï 1997; Shepel *et al*. 1998; Sen and Churchill 2001; Cui *et al*. 2006), aimed to map “count” QTL that determine a unique aspect of biological systems.

There is no difficulty in using current-interval mapping, composite-interval mapping, or multiple-interval mapping approaches to map the phenotypic plasticity of a biological trait, defined as the phenotypic difference of the same genotype across a range of environments (Wu 1998; Leips and Mackay 2000; Kliebenstein *et al*. 2002; Ungerer *et al*. 2003), although a mechanistic basis of the difference of QTL expression between different environments is not considered. However, by extending these mapping approaches into a multivariate case in which phenotypic vales of a trait in different environments can be regarded as different “traits” (Jansen *et al*. 1995; Jiang and Zeng 1995), environment-dependent genetic effects of a QTL can be estimated and tested. The idea of multitrait QTL mapping was used in this study to map a QTL that triggers a pleiotropic effect on a count trait expressed in different environments. If such a pleiotropic effect is different between environments, genotype-by-environment interactions at this QTL result (Ungerer *et al*. 2003).

Our model capitalizes on the merit of a randomized complete block design, characteristic of the study of phenotypic plasticity, in which the same genotypes are grown in different environments with multiple replicates per environment used to minimize microenvironmental noises. By incorporating a multivariate Poisson distribution, our model allows the modeling of count data in individual replicates, thus avoiding losing the integer-valued structure of the data when means over replicate traits are taken. The estimation of the multivariate Poisson distribution is a difficult issue, but a powerful EM algorithm has been implemented to estimate the parameters (Karlis 2003; Tsiamyrtzis and Karlis 2005). We integrate this multivariate Poisson-based model into a standard mixture model for QTL mapping, leading to a two-stage hierarchical EM algorithm that can estimate and test the genetic effects of a QTL expressed in different environments. Through computer simulation, we investigated the statistical behavior of the new model in terms of its estimation accuracy and precision and the power to detect environment-dependent QTL for a count trait. Our simulation results suggest that a modest number of mapping progeny (100) can offer reasonable parameter estimation and power that is improved with increasing genotypic differences. The application of our model to a real example for Populus genetics leads to the detection of a significant QTL for the number of sylleptic branches that increases significantly with better fertilization. By the continuous development of an axillary meristem into a branch without an intervening stage of rest, syllepsis has been thought to be an important determinant of wood production (Wu and Hinckley 2001). The identification of an environment-sensitive QTL for the number of sylleptic branches is consistent with previous findings that this type of branch is under strong genetic control and highly plastic to changes in environment (Wu and Stettler 1998).

Our model can be used to test one of three existing hypotheses proposed to explain phenotypic plasticity—allelic sensitivity; *i.e*., different expression of a QTL in different environments causes the phenotypic plasticity of a trait (Via and Lande 1985). Although the pseudotest backcross design used, in which only two, rather than all three possible, genotypes exist, does not allow for the test of the overdominance hypothesis (Gillespie and Turelli 1989), this hypothesis can be tested by including F_{2}-type codominant markers and QTL, *i.e*., those loci that are segregating for both parents and thus form two homozygotes and one heterozygote, in the mapping population derived from outcrossing trees (Lin *et al*. 2003). It is also possible that our model can be modified to model the gene regulation hypothesis of phenotypic plasticity by assuming different but epistatically interacting QTL for plastic responses and means of a trait across different environments (Scheiner and Lyman 1989; Weber and Scheiner 1992). In addition, there are some other issues that deserve further investigation. In genetics, it is practically useful to implement composite-interval mapping (Zeng 1994; Jansen and Stam 1994) and multiple-interval mapping (Kao *et al*. 1999) into our Poisson-based mapping model because these techniques have a capacity to remove the background noise of QTL effects and/or map a network of interacting QTL at the same time. In the current statistical model, we used a single Poisson parameter (θ_{0}) to model the covariance structure among different replicates, but the power of our model can increase if more sophisticated modeling is used for the covariance structure (Karlis 2003). Also, Cui *et al*. (2006) recently generalized the Poisson-based mapping model to take into account the dispersion of count data. By estimating and testing a dispersion parameter in our model, we are able to estimate the genetic effects of a QTL on the direction and degree of the dispersion of a count trait.

Characterizing how the genetic architecture of a complex trait differs across environments is an important first step toward elucidating the mechanistic basis of genotype-by-environment interactions and ultimately predicting final phenotypes in a given environment. To achieve this goal, it is necessary to identify genes underlying quantitative variation in plastic response of various complex traits (including count traits) to changes in environments and determine how these genes act singly or interact epistatically to determine different phenotypes across ecologically relevant environments. The model proposed in this article and its extensions as discussed above will provide a powerful tool to gain insights into how genomes and environments interact to determine phenotypes and further alter the evolutionary process of adaption.

## Acknowledgments

We thank the graduate students in STA 6178 class of Spring 2007 at the University of Florida for their contributions to this work. The preparation of this manuscript was supported by National Science Foundation grant 0540745 to R.W.

## Footnotes

↵1 These authors contributed equally to this work.

Communicating editor: G. Gibson

- Received September 11, 2007.
- Accepted February 14, 2008.

- Copyright © 2008 by the Genetics Society of America