Genetic Regulation of Phenotypic Plasticity and Canalisation in Yeast Growth

The ability of a genotype to show diverse phenotypes in different environments is called phenotypic plasticity. Phenotypic plasticity helps populations to evade extinctions in novel environments, facilitates adaptation and fuels evolution. However, most studies focus on understanding the genetic basis of phenotypic regulation in specific environments. As a result, while it’s evolutionary relevance is well established, genetic mechanisms regulating phenotypic plasticity and their overlap with the environment specific regulators is not well understood. Saccharomyces cerevisiae is highly sensitive to the environment, which acts as not just external stimulus but also as signalling cue for this unicellular, sessile organism. We used a previously published dataset of a biparental yeast population grown in 34 diverse environments and mapped genetic loci regulating variation in phenotypic plasticity, plasticity QTL, and compared them with environment-specific QTL. Plasticity QTL is one whose one allele exhibits high plasticity whereas the other shows a relatively canalised behaviour. We mapped phenotypic plasticity using two parameters–environmental variance, an environmental order-independent parameter and reaction norm (slope), an environmental order-dependent parameter. Our results show a partial overlap between pleiotropic QTL and plasticity QTL such that while some plasticity QTL are also pleiotropic, others have a significant effect on phenotypic plasticity without being significant in any environment independently. Furthermore, while some plasticity QTL are revealed only in specific environmental orders, we identify large effect plasticity QTL, which are order-independent such that whatever the order of the environments, one allele is always plastic and the other is canalised. Finally, we show that the environments can be divided into two categories based on the phenotypic diversity of the population within them and the two categories have differential regulators of phenotypic plasticity. Our results highlight the importance of identifying genetic regulators of phenotypic plasticity to comprehensively understand the genotype-phenotype map.


Introduction
A single genotype cannot have high fitness in all conditions. Instead different genotypes show varying degrees of fitness in different environments, and therefore phenotype of a genotype is dependent on the environment. The ability of a single genotype to show different phenotypes in different environments is called phenotypic plasticity [1]. On the other hand, ability of a genotype to show the same phenotype independent of the environment is termed as canalisation [2]. Phenotypic plasticity facilitates adaptation to novel environments by allowing the population to exhibit a diverse range of phenotypes [3]. It is ubiquitous in nature and shown to be a major force in adaptation, be it adaptation to climate change, altitude, nutrition, multi-cellularity, etc. [4,5]. Consequently, phenotypic plasticity is one of the major drivers of evolution [2,6,7].
During adaptation, stabilising selection acts on the population such that the phenotype gets stabilised or canalised within an environment and across multiple environments [8]. One of the ways this canalisation is proposed to get perturbed is when this adapted population encounters a novel or rare environment. This perturbation of canalisation allows the population to exhibit a range of phenotypes thus facilitating adaptation. Canalisation and plasticity are dynamic, mutually dependent processes and a population switches between these two states depending on the environments encountered [9,10]. While a canalised phenotype would be beneficial in environments to which the population has adapted to, a plastic phenotype would be advantageous in a novel or rare environment [6]. Hence the same genotype is capable of showing a canalised or plastic behaviour depending on the environments considered and different genetic regulators may regulate phenotypic plasticity in varying environments.
While the importance of plasticity in adaptation and evolution has been established by multiple studies [2,11,12], these studies are mostly either theoretical or conducted in naturally occurring populations. Therefore, while evidence for phenotypic plasticity has been documented in multiple organisms across diverse phenotypes, its genetic regulation is not clearly understood. Additionally, most studies that attempt to understand the genetic regulation of a phenotype focus on either a single environment or multiple environments independently [13][14][15]. As a result, while our knowledge about genetic regulation of a phenotype in different environments is fairly comprehensive, we do not understand the genetic regulation of plasticity and canalisation across diverse environments. While phenotypic plasticity is mainly invoked to study the adaptability of natural populations, its ubiquity and role in evolution indicates that it should also be important for understanding the genetic architecture of complex traits [16].
Quantitative trait locus (QTL) mapping provides a good way to identify regulators of phenotypic plasticity. Phenotypes of most loci show environment dependence [17]. By this definition, all loci showing gene-environment interaction (GEI) exhibit phenotypic plasticity. However, a plasticity QTL is a locus whose one allele shows a canalised behaviour whereas the other allele shows phenotypic plasticity across diverse environments [18] (Fig 1A and 1B). If two genetically diverse strains have encountered and adapted to varied environments, or adapted to the same environments using different mechanisms, then crossing these strains will disrupt these mechanisms and allow identification of loci with differential plasticity in this biparental population.
While multiple studies have performed QTL mapping to identify plasticity QTL, they were either done across pairs of environments or continuums of environments [18,19]. However, in nature, populations encounter diverse environments, capable of affecting the phenotype, either simultaneously or consecutively. In-lab evolution studies have shown that the order of environments encountered during the course of evolution can dictate which alleles eventually get fixed in a population [20]. Parallel to this, it is probable that the order of encountering these environments would determine the plasticity of the genotype, which would in turn determine the selection forces that act on it (Fig 1). Different genotypes can show different ranges of phenotypic plasticity depending on the order of the environments and different parameters are required to capture the plasticity in different environmental groups (Fig 1C and 1D). Hence in order to comprehensively identify the regulators of phenotypic plasticity, both diversity of environments and their order should be considered.
In this paper, we asked the following questions: can we identify plasticity QTL across a large number of heterogeneous environments? How do these plasticity QTL respond to different types and orders of environments? Finally, what is the association between pleiotropic regulators of the phenotype and plasticity regulators? Are loci that regulate plasticity and that are pleiotropic across multiple environments same such that all pleiotropic loci contribute to plasticity, or are these loci different and hence not identified in environment-specific QTL mapping?
S. cerevisiae provides an ideal system to identify the genetic regulators of phenotypic plasticity, since environment serves as both external stimulus as well as signalling cue for this unicellular, sessile, organism. Yeast growth is highly responsive to environments and has been shown to be differentially regulated in different environments [17,21,22]. In this study, using growth phenotype measured in 34 diverse environments for a large yeast biparental population [14], we measured phenotypic plasticity using two statistics: an environmental order-independent statistic-Environmental variance (Var E ), and an environmental order-dependent statistic, Sum Genotype A1 and A2 are represented in blue and red colours respectively. Var E refers to environmental variance whereas X Slope refers to sum of slopes, as described in Methods. y-axis denotes the phenotype and x-axis denotes discrete environments arranged in different orders. (A) Genotype A1 and A2 have significant differences in multiple environments but are both equally plastic. (B) A1 is plastic and A2 is canalised. (C) and (D) shows the same environments arranged in different orders which have no effect on environmental variance but have different impact on reaction norms or sum of slopes. of slopes (reaction norms) ( X Slope) (Fig 1). Fig 1 shows that both these parameters capture different aspects of phenotypic plasticity. Fig 1A shows that genotypes with difference in phenotype across diverse environments do not necessarily have differential plasticity; 1B shows two genotypes with differential plasticity; and Fig 1C and 1D show that while the environmental order has no bearing on environmental variance, the value of the reactions norms is highly sensitive to the order of the environments encountered. We used these two parameters to identify loci with differential effects on phenotypic plasticity, plasticity QTL. To the best of our knowledge, this is the first study to identify genetic regulation of phenotypic plasticity and canalisation across such a diverse set of environments. These genetic regulators of phenotypic plasticity may play an important role in explaining missing heritability and understanding the genetic regulation of complex traits especially human disease that are influenced by multiple environmental conditions.

Dataset
The raw growth data used in this study was derived from a previously published study by Bloom et al. [14], in which the experimental procedures are described in detail. The data we used was generated for 1,008 segregants derived from a cross between yeast strains BY (a laboratory strain) and RM11-1a (a wine isolate, indicated as RM). These segregants were genotyped for a total of 11,623 polymorphic markers and were grown and phenotyped for colony size in 46 different conditions. Of these 46 conditions, we selected 34 conditions based on following three criteria: (i) segregant phenotype should show normal distribution; (ii) environments should be closer or mimic naturally occurring environmental conditions; (iii) since different degrees of the same environmental stresses can invoke correlated phenotypes, biasing our analysis, the environments should be heterogeneous and not continuums. This filtering removed environments like rapamycin, pH and temperature gradients, etc.

Single QTL Mapping
QTL mapping was carried out as described previously [22]. In brief, the R/qtl package [23,24] was used to identify QTL separately for colony size in each environment. QTL were identified using the LOD score, which is the log 10 of the ratio of the likelihood of the experimental hypothesis to the likelihood of the null hypothesis [24]. An interval mapping method ('scanone' function in R/qtl) was used to compute this LOD score using the Haley-Knott regression algorithm [23].
The following formula was used to calculate the F-score, which was further used to derive the LOD score. At a particular marker, let segregant i's phenotypic value be y ij where j can take two values (j = 1: BY allele and j = 2: RM allele).

ðN À kÞ
here, N is the total number of segregants, n 1 and n 2 are the number of segregants having the BY and RM allele respectively (k = 2) and y i is the genotypic mean of allele j.
Let df denote the degrees of freedom (df = 1 for a backcross and df = 2 for an intercross). The LOD score is accordingly derived as follows: Under the null hypothesis, there is no significant difference in the means at the marker under consideration while under the alternative hypothesis, there is a presence of a QTL.

Plasticity QTL Mapping
Plasticity QTL mapping was performed using the same methodology as described for QTL mapping, using environmental variance and sum of slopes as phenotypes, instead of colony size.
Environmental variance (Var E ) was computed for each segregant separately for high (Hv) and low (Lv) variance environments: where, x is phenotype of a segregant in an environment, μ is the average phenotype across n environments. n = 10 for Hv and n = 24 for Lv environments. For mapping in sub-groups of Hv environments, n was 3 and 4, respectively.
Sum of slopes ( X Slope) was calculated for each segregant for each order of environments using the following formula: Where n is number of environments in a given order, x is the phenotype in the environment and c is the constant that represents difference between the two environments. Since all the environments are heterogeneous discrete environments and do not represent a continuum, the difference between them is always a constant, thus c was given a value of 1.

Random orders and allele specific plasticity QTL
Environmental order for calculating the sum of slopes was determined in two ways for both Hv (10 environments) and Lv (24 environments) environments: random orders and allele specific orders. In random orders, 10 random orders of environments were generated. For a particular order, each environment was given a single unique position, such that there were no repetitions of environments. For each order, sum of slopes was calculated for all segregants and QTL mapping was performed. In allele specific orders, orders of environments were generated separately for both BY and RM alleles for each marker. For each allele at a particular marker, the environments were ordered such that the mean of the segregants carrying that allele have the least possible sum of slopes. In other words, the mean of the population is canalised across the environmental order. Sum of slopes was calculated for this order for all segregants and QTL mapping was performed. Therefore QTL mapping was performed 6 times using sum of slopes for different types of environmental orders = 2 (Hv and Lv random orders) + 2 (Hv and Lv allele specific order with BY allele canalised) + 2 (Hv and Lv allele specific orders with RM allele canalised). Total number of environmental orders tested for each group = 10 (random orders) +~11,623 (BY allele canalised at each marker) +~11,623 (RM allele canalised at each marker).

Environments fall into two categories based on the variance of the segregants
In the previously published dataset [14], we computed the variance of all segregants across 34 environments to identify the range of phenotypic plasticity exhibited by the individuals of the population. A high variance would indicate high diversity of the phenotype of the segregant across the environments (high phenotypic plasticity) whereas a low variance would suggest similar phenotype across all environments (canalisation). The phenotypic variance showed a normal distribution indicating that it was a complex trait with a fraction of individuals showing highly canalised and highly plastic behaviour (Fig 2A and S1A and S1B Fig). There was no association between the variance and average phenotype of the segregants (R 2 = 0.0007) indicating that segregants with both high and low average phenotype could show high variance. Apart from the genotype, the environments considered also determine the plasticity of an individual. We have previously shown that while a population shows highly buffered phenotype in one environment, this buffering can be lost in others [25]. Hence, we compared the phenotypic variance of the segregants within each environment ( Fig 2B). The variance in the 34 environments did not show either a normal or a bimodal distribution but a highly left skewed distribution with a median of 4.2 ( Fig 2B). Hence we categorised the environments that were within the first quartile (0 to 8) in the category Lv environments. While the remaining 10 environments showed a large range of variance, splitting them into smaller number of environments could have reduced the statistical significance of the variance and slope phenotypes. Therefore, we categorised these 10 environments as Hv environments (Fig 2B). We calculated variance of each segregant in Lv and Hv environments independently, and found no correlation between the two values ( Fig 2C). This indicates that a segregant with highly variable phenotype in Lv environments can be either plastic or canalised in the Hv environments and vice versa. We also calculated the mean of each segregant across Hv and Lv environments, and found it to be poorly correlated (R 2 = 0.03, S2A Fig). Furthermore, if genetic regulation between random sets of Lv environments was as diverse as that between Hv and Lv environments, then we should observe poor correlation among Lv environments. We sampled two random sets of 10 environments each from the Lv category and computed correlation of mean values of segregants. These two sets had non-overlapping environments such that the presence of common environments does not bias the correlation. We observed a significantly high correlation between mean across these two sets (R 2 = 0.38, P < 0.01, S2B Fig), which indicated similar genetic regulation in Lv environments, but differential regulation across the Hv and Lv environments.

Different loci are pleiotropic in high and low variance environments
Studies have shown that while most yeast growth QTL tends to be environment specific, some loci have pleiotropic effects. A pleiotropic locus is one that has an effect on the phenotype across multiple environments. In order to determine whether plasticity QTL are the same as, or a subset of, or entirely different from pleiotropic QTL, we carried out QTL mapping in each environment (see Methods). A complete overlap of the large effect QTL and a high overlap of small effect QTL was observed between this study and the original study by Bloom et al. [14] (S1 Table) reconfirming our mapping results. We first compared the pleiotropic loci identified in multiple environments. A locus was designated as pleiotropic if it has an effect in 4 or more environments with a LOD peak within 40kb interval across these environments. Multiple QTL were identified to be pleiotropic across the 34 environments (Table 1).
We next asked if the pleiotropic loci were different between the Hv and Lv environments. We found that some pleiotropic loci were common, but others were specific to only Hv or Lv environments (Fisher's Exact test P < 0.1, Table 1 and S1 Table). This shows that there exists a difference in genetic regulation of the phenotype between the Hv and Lv environments, as predicted by poor correlation of mean across Hv and Lv environments but strong correlation among Lv environments (S2 Fig). Previously published fine mapping studies done using the BYxRM segregant populations provide potential candidate genes in many of these loci. chrXIVb and chrXVa peaks have been identified in multiple environments and fine-mapped to pleiotropic genes like MKT1 [13] and IRA2 [13,26] respectively, however in this study neither of these were identified as plasticity QTL in either category of environments. Another pleiotropic QTL, chrXIII locus has been previously associated with yeast chronological lifespan and  telomere length with gene BUL2 as causative [27]. Finally, chrV QTL effected colony morphology with GPA2 as causal gene [28]. While chrXIVa QTL has not been fine-mapped to any gene, various peaks identified in single QTL and plasticity QTL mapping (see below) indicated that causal gene could be KRE33, a protein required for biogenesis of small ribosomal subunit with its human homolog implicated in several types of cancer and premature ageing [29]. However, confirmation of involvement of these candidate causative genes is pending experimental validation.

Identifying plasticity QTL using environmental variance
In order to identify plasticity QTL, the first step is to determine a parameter that captures plasticity of segregants. We used modifications of two commonly used parameters: variance and reaction norm or slope [18,30]. Commonly applied data normalisation across environments enhances the power of comparing effect of loci across two environments and helps identifying GEI. However, it also makes the allelic effects symmetric thereby making both alleles equally plastic which results in an inability to distinguish between plastic and canalised alleles ( Fig 1A). Therefore, since the aim of this paper was to identify plasticity QTL and not GEI, we normalised the phenotype within an environment but not across environments. While this reduced the power of identifying QTL, the ability to identify plasticity QTL was preserved. Whether one does across-environment normalisation or not, this has no bearing on the QTL identified within an environment [17]. Environmental variance (Var E ) refers to the variance of the phenotype of a segregant across multiple environments. As discussed above, high variance would indicate that the segregant has diverse or plastic phenotype across environments and low variance would suggest that the segregant shows similar phenotype, or canalised behaviour, across environments. Since the scale of variance was different for Hv and Lv environments (Fig 2B), Var E was calculated for each segregant independently for each class of environments. As a result, we got two phenotypes for each segregant: Var E in Hv and Var E in Lv environments. We observed no correlation between average phenotype and segregant Var E , indicating that the two properties were not significantly related (Pearson correlation P > 0.1). We then performed QTL mapping for these two phenotypes. While the overall LOD scores identified were lower than conventional single environment QTL mapping, the peaks were significant (Fig 3A and 3D and S2 Table, permutation P < 0.01). Two peaks were identified in Hv (Fig 3B and 3C) and one in Lv environments ( Fig 3E) with a LOD score > 2.0 (P < 0.01). The highest peak in Lv environments, chrXIVa locus was pleiotropic and unique to this class of environments (Table 1). One peak in Hv environments was pleiotropic (chrXIII locus) whereas the other was not (chrV locus). Interestingly, for both the peaks in Hv environments, on chrV and chrXIII, the RM allele had higher environmental variance than BY allele, whereas for the peak in chrXIVa locus in Lv environments, the BY allele showed higher environmental variance (Fig 3 and S2 Table). Surprisingly in single QTL mapping, BY allele of chrXIVa that is a more plastic allele, had lower mean than the RM allele in almost all cases.
While the environments with variance greater than 8 were categorised as Hv environments, as the Fig 2B shows, the highest variable environments show large variance values and can possibly themselves be split further into two subgroups. Therefore, we split 7 Hv environments (variance greater than 20) into two subgroups-Hv_subgroup1, consisting of 4 environments, and Hv_subgroup2, consisting of 3 environments (see S1 Table). Var E was calculated for each segregant independently for each subgroups and QTL mapping was performed as previously discussed (S3 Fig and S2 Table). Some peaks were specific to each subgroup, for example, a peak on chrX was specific to Hv_subgroup1 and another on chrXII was specific to  Hv_subgroup2 (S3 Fig). However, the large effect chrXIII locus that was both pleiotropic and plastic in all Hv environments, was also identified in both the subgroups (S3 Fig and S2 Table) supporting the original categorisation of Hv and Lv environments Many loci that were pleiotropic across different environments were not identified as plasticity QTL. A stark example is the chrXIVb locus that has been identified as a pleiotropic locus in many environments but had no effect on phenotypic plasticity (Table 1).

Identifying plasticity QTL using sum of slopes
While Var E provides an unbiased measure of phenotypic plasticity, it is not sensitive to relatively small changes in the phenotype (Fig 1D). As a result, most GEI studies calculate reaction norms or slopes to identify small effect but significant changes in the phenotype across environments. Usually GEI analysis is performed for a pair of environments [17,22]. As shown by these studies, the steeper the slope of the reaction norm, the more plastic is the genotype. While sensitive, this method can be used only for 4-5 environments or continuums of environments. Large number of heterogeneous environments results in multiple pairwise comparisons that are difficult to both compute and compare. We overcame this shortcoming by computing a novel parameter called sum of slopes ( X Slope, see Methods, Fig 1). Briefly, we arrange the environments in different orders and calculate slopes between consecutive environments. So that slopes in opposite direction do not cancel each other, absolute values of these slopes are summed to obtain a value of the parameter. Higher the sum of slope value, more plastic is the individual. Unlike Var E , sum of slopes will depend upon the order of the environments considered (Fig 1C and 1D). We asked the following questions: how much overlap will be observed in the plasticity QTL mapped using these two different parameters? Will identification of plasticity QTL using sum of slopes depend on the order of the environments?
As done for Var E , we calculated sum of slopes for each segregant separately for the Hv and Lv environments. For each category, we used two different strategies to compute the order of the environments. First strategy was to generate random orders where, using permutations, we computed 10 random orders of the environments and then calculated sum of slopes for each segregant for an order and used this as a phenotype for mapping. As a result, we obtained plasticity QTL for each order of the environments, for both Hv and Lv environments separately (S3 Table, permutation P < 0.01). Second strategy was to generate allele specific environmental orders, which takes into consideration that different alleles might have evolved as a result of different selection pressures and hence show canalisation across different orders of environments. While 10 combinations is a substantial number, it may not be exhaustive enough to identify canalisation orders for all alleles. Therefore, we ordered the environments for each allele of each marker independently. For both Hv and Lv environments independently, for each locus, the environments were ordered to have the least possible sum of slopes for one allele. This order was then used to calculate sum of slopes for all segregants and the values were used for plasticity QTL mapping for that particular marker. The same was done for the other allele separately. This was done for all markers independently. Therefore, the total number of environmental orders tested was equal to the product of number of markers, two categories of environment and two alleles. Thus, the QTL were mapped for a canalised mean of each allele for each marker, in both categories of environments (Fig 4 and S4 Table, permutation P < 0.01).
Higher LOD scores and larger number of plasticity QTL were identified for sum of slopes than environmental variance (Table 1 and S2 and S3 Tables). For random order analyses, the plasticity QTL identified depended on the order of the environments. We compiled the results to identify peaks that were identified in most environmental orders. Certain plasticity QTL were identified in more than half of 10 random environmental orders, i.e. they were independent of the environmental order (Table 1). While 4 peaks were identified in majority of the environmental orders consisting of Lv environments, only a single peak was consistently identified in Hv environments (Table 1 and S3 Table). These loci included the ones identified using Var E , as well as unique to sum of slopes (Table 1). Distinct sets of peaks were identified in Hv and Lv environments using allele specific environmental orders (S4 Table). Moreover, LOD scores were higher for Lv than Hv probably due more noise in the phenotype, sum of slopes, in Hv compared to Lv environments (S3 and S4 Tables). Additionally, like the plasticity QTL identified depended on the random order, the identification of the plasticity QTL using allele specific order depended on the allele whose mean was canalised (Fig 4 and S4 Table). However, we also identified plasticity QTL that were independent of the allele whose mean effect was canalised, i.e. they were identified independent of whether the RM or BY allele was canalised. These overlapped with the plasticity QTL that were identified in most random orders of environments (Table 1).
We compared plasticity QTL identified using three strategies: Var E , sum of slopes with random orders and sum of slopes with allele specific orders ( Table 1). As proposed in Fig 1, both Var E and sum of slopes are capable of identifying differences in plasticity to different extents and measuring both of them is required to identify the genetic regulators of phenotypic plasticity. While several QTL were specific to the parameter or environmental order used, two loci chrXIII in Hv and chrXIVa in Lv environments were identified in all three methods ( Table 1). Identification of these plasticity QTL through independent strategies emphasises their definite ability to regulate phenotypic plasticity.
Comparison of sum of slopes revealed that, as expected, the value of this parameter was less for Lv than for Hv environments. However, canalisation of mean of the allele, i.e. the lowest sum of slopes of mean, as done for allele specific order, did not necessarily result in reduced sum of slopes of the segregants carrying the allele (S4 Table). For plasticity QTL that were identified independent of the allele, the same allele had higher sum of slopes of segregants independent of the allele whose mean was canalised (S3 and S4 Tables). This explains why some plasticity QTL were identified irrespective of the environmental order. Furthermore, this shows that canalisation of the population mean does not always reflect canalisation of the individuals within the population (Fig 5A and 5B). An allele can have a canalised mean but differential plasticity of individuals. This was observed for both Hv and Lv environments (Fig 5A  and 5B). Furthermore, our results show that while environmental order can uncover the difference in plasticity between two alleles, a canalised allele will always be canalised independent of the environmental order (S4 Table).
High variance of sum of slopes within an allele would indicate diversity of phenotypic plasticity. While there was no association between mean and variance of segregant values across environments, we found that there was a positive association between the mean and variance of sum of slopes between various alleles in both Hv and Lv environments indicating that the allele with higher sum of slopes also showed more diversity (Fig 6). Therefore, the segregants carrying the more plastic allele did not show same pattern of phenotypic plasticity but demonstrated a diversity of patterns, potentially to facilitate adaptation to diverse environments. Hence, our results show that the more plastic allele also results in revelation of hidden reaction norms.

Discussion
Our study identifies loci with differential effects on phenotypic plasticity in heterogeneous environments. We show that regulation of phenotypic plasticity is overlapping but different than the regulation of phenotypic variation in each environment. This has implications not only on adaptation and evolution, but also on understanding the genetic architecture of genotype-phenotype map. While different plasticity QTL were identified using different parameters of plasticity and in different environmental orders, some of these plasticity QTL were identified in all mapping methods indicating their robust role in regulating phenotypic plasticity.
Phenotypic plasticity is a property of the genotype, unveiled by the environments. We show that environments can be divided into two categories based on phenotypic variance of the population termed here as Hv and Lv environments (Fig 2). Such a distinction has been hypothesised by previous studies [4], which propose that when a population is adapted to a particular environment, then stabilising selection acts on the population, such that most individuals of the population show similar phenotype which is close to the fitness optimum (low variance). When the population encounters a novel or rare environment, this buffering is perturbed releasing high diversity of individual phenotypes (high variance), which can facilitate adaptation. In the light of this current evolutionary understanding of plasticity and canalisation, we infer our results from a biparental population as follows: the Lv environments are the ones in which either one or both strains have adapted to in the course of their evolutionary history whereas the Hv environments are potentially novel environments [31]. This conclusion is further facilitated by identification of different QTL as well plasticity QTL in both these categories of environments (Table 1). Differential enrichment of pleiotropic QTL in the two categories indicates a common regulation of the phenotype in the canalised or Lv environments. Additionally, disruption of canalisation in the recombinant population may explain why the large effect and consistent plasticity QTL were identified in Lv than the Hv environments. Genetic recombination disrupts the evolved canalisation mechanisms therefore resulting in identification of plasticity QTL in Lv environments, whereas poor or no canalisation mechanisms exist for Hv environments, which results in high plasticity of all alleles. This results in reduced LOD score of plasticity QTL identified.
As proposed in Fig 1, our results show that plasticity QTL are not same as pleiotropic QTL. Almost all loci show GEI and large effect pleiotropic loci show large effect GEI [17]. However, we observed only a partial overlap between pleiotropic QTL and plasticity QTL. While some large effect QTL (like chrXIII and chrXIVa) also had pleiotropic effects, others like chrV and chrXIVc did not show pleiotropy but were equally significant plasticity QTL. In fact, while chrXIVa and chrXIVc were in a relative close physical distance, within 160kb (Table 1), they had opposite effects on plasticity of the alleles: BY allele of chrXIVa showed high plasticity and RM allele of chrXIVc showed high plasticity (S3 and S4 Tables). This indicates that genetic regulation of phenotypic plasticity is overlapping, but different than genetic regulation within each environment. This further emphasises that in order to understand the genotype-phenotype map and the function of identified molecular regulatory hubs, it is important to not only understand their effects in one environment or phenotypes but across different environments.
In a previous study, we showed the biological implication of mean and variance of a population [25]. We showed that higher variance was associated with phenotypic manifestation of cryptic or hidden variants. Additionally, a high phenotypic variance could either be associated with a higher or a lower phenotypic mean depending on the environment. Here we show a strong correlation between mean and variance of phenotypic plasticity (Fig 6). Interestingly, in both Hv and Lv environments, the allele with a higher mean of plasticity also had a higher variance (Fig 6). This indicates that segregants containing the more plastic alleles exhibit a diverse range of phenotypic plasticity, potentially to facilitate adaptation in diverse environmental conditions. The high variance of plasticity values (both Var E and sum of slopes) suggests presence of epistasis, resulting in revelation of hidden reaction norms [4] or cryptic genetic variants with diverse effects across environments. Along with shedding light on mechanisms of regulation of phenotypic plasticity, this suggests an association between genetic regulation of cryptic genetic variation and phenotypic plasticity [32].
In conclusion, by identifying genetic regulators of phenotypic plasticity and canalisation, our results highlight that genetic regulation of a phenotype in an environment may depend not only upon mechanisms directly evolved in that environment but maybe a result of evolution in a diverse range of environments [16,33]. While commenting on the evolutionary nature of the identified plasticity QTL is beyond the scope of our results, our study opens new avenues of exploring population genetic data and understanding the underlying basis of the genetic architecture. Differential regulation of phenotypic plasticity provides a potential reason underlying the high interconnectivity observed in the genotype-phenotype map. This interconnectivity could be an outcome of cross talk between different genetic modules that either maintain canalisation or induce plasticity across different environments and phenotypes. This has profound implications, especially on understanding adaptation mechanisms in naturally occurring plant and animal populations, development [34] as well as understanding the molecular basis of regulation of complex human diseases highly susceptible to environmental conditions [35] such as metabolic and psychological disorders. x-axis shows the difference between mean value of sum of slopes of alleles for different peaks, BY-RM, and y-axis refers to difference between variance of sum of slopes of alleles, BY-RM. See S3 Table for