Large-effect flowering time mutations reveal conditionally adaptive paths through fitness landscapes in Arabidopsis thaliana.

Contrary to previous assumptions that most mutations are deleterious, there is increasing evidence for persistence of large-effect mutations in natural populations. A possible explanation for these observations is that mutant phenotypes and fitness may depend upon the specific environmental conditions to which a mutant is exposed. Here, we tested this hypothesis by growing large-effect flowering time mutants of Arabidopsis thaliana in multiple field sites and seasons to quantify their fitness effects in realistic natural conditions. By constructing environment-specific fitness landscapes based on flowering time and branching architecture, we observed that a subset of mutations increased fitness, but only in specific environments. These mutations increased fitness via different paths: through shifting flowering time, branching, or both. Branching was under stronger selection, but flowering time was more genetically variable, pointing to the importance of indirect selection on mutations through their pleiotropic effects on multiple phenotypes. Finally, mutations in hub genes with greater connectedness in their regulatory networks had greater effects on both phenotypes and fitness. Together, these findings indicate that large-effect mutations may persist in populations because they influence traits that are adaptive only under specific environmental conditions. Understanding their evolutionary dynamics therefore requires measuring their effects in multiple natural environments.

Throughout much of the 20th century, mutations were assumed to be largely deleterious, diminishing the adaptive value of a trait in a particular environment (1–3). However, more recently, a robust empirical body of work has shown that the frequency of beneficial mutations considerably exceeds previous assumptions (4–10). Many mutations, even of large effect, are known to have become resident in populations and persist as segregating alleles (11–16). The few experiments testing mutations’ fitness effects in nature have revealed complex patterns among environments, with effects in one environment poorly predicting them in another (6, 17–19). A major reason for this complexity may be that environments affect trait expression, but few, if any, studies have gauged the mutation–phenotype relationship in multiple natural environments. Rather than examining unbiased mutation spectra consisting of many neutral mutations present in mutation accumulation lines (19), one way to test this is to examine targeted mutations known to affect ecologically important phenotypes. Here, we do this by constructing genotype–phenotype–fitness landscapes of Arabidopsis thaliana mutants harboring large-effect mutations targeted to known flowering-time genes in common gardens across its European climate range.

Fitness landscapes were introduced as a metaphor to conceptualize how populations evolve through states of low fitness to high fitness (20) and have since been extended and formalized to express mathematically the relationships among genotypes, traits, and fitness (21–23). In practice, empirical fitness landscapes often examine only one-half of the genotype–phenotype–fitness relationship. In the phenotype–fitness approach, landscapes relate trait values to fitness so that a peak represents a combination of trait values that is highly fit. However, these landscapes are usually constructed without reference to specific genes that underlie these traits (24–36). Genotype–fitness landscapes, on the other hand, relate combinations of alleles to fitness so that a peak represents a highly fit genotype. While powerful for understanding evolution, measurements of genotype–fitness landscapes have been rare in multicellular organisms, and most do not consider higher order traits directly related to fitness (20, 37–50).

A major goal in evolutionary biology is to synthesize these approaches in a complete genotype–phenotype–fitness landscape. However, in eukaryotic organisms, combining phenotype–fitness and genotype–fitness landscapes is frequently limited by either a lack of knowledge of adaptive traits or of genes that control them. The genetic model species A. thaliana provides a valuable system to tackle this challenge: The genetic control of ecologically important phenotypic traits is well studied (51–54), and field experiments allow the measurement of phenotypic selection in natural environments (55–61).

Both flowering time and branching architecture are potentially important adaptive traits in A. thaliana and other annual plants. The timing of flowering determines the critical transition from vegetative to reproductive growth; aligning this transition with favorable environmental conditions is critical to successful reproduction. Flowering time varies clinally across the species’ range (53, 56, 60, 62–66), evidence for adaptive differentiation in this trait across different environments. This clinal differentiation is due to natural variation in responses to environmental signals like temperature and day length (66–69). Flowering time is closely linked to branching architecture through development. Prior to flowering, Arabidopsis plants produce rosette leaves and rarely secondary branches. However, after commitment to flowering, side branches are produced from axillary meristems, and these side branches also produce flowers (70, 71). The number of axillary meristems and branches available to produce seeds depends upon flowering time since prolonged vegetative growth due to delayed flowering allows a plant to develop more leaves and axillary meristems (72–74). Thus fecundity is at least partially a function of flowering time-dependent branching architecture (75). Importantly, the hormone gibberellin (GA) partially decouples branching from fecundity since it promotes branching but simultaneously inhibits flower formation (76). This decoupling enables branching and reproductive phenology to act as partially independent developmental modules that may be independently selected upon but that are still partially united. Thus, plant fitness depends on a balance between branching and the timing of flowering to determine fitness in morphological and phenological space. Since both flowering time and gibberellin flux depend upon growing conditions (77–79), this balance is likely to shift among environments.

The genotype–phenotype map is further complicated by the complex and nonlinear relationships among genes, especially for quantitative traits like branching and flowering time with highly polygenic bases. Genes influencing these traits interact through multilevel networks. Some are hubs that influence a large number of other genes, but others are subject to multiple inputs. These highly connected genes within regulatory networks may allow mutations in them to be buffered against producing pernicious phenotypes (80–82). The degree of this buffering is a function of both a gene’s network connectedness and redundancy permitting increases or decreases in these factors, respectively, to lead to more severe deviation in phenotype when mutated (80, 83–86). Together, this implies that potentially adaptive mutations may be masked due to network topography, robustness, and redundancy. Thus, adaptive zones in fitness landscapes may be accessible only through mutations in large-effect hub genes. On the other hand, such robustness is also a function of the environment so that mutations buffered in one environment may be exposed in another (87–90). In this way, genetic network position can interact with environment to shift mutants in fitness landscapes.

Although the accumulation of many small-effect mutations can collectively alter phenotypes, we focus on individual large-effect allelic variants induced by chemical and radiation mutagenesis. Single, large-effect mutations represent a crucial evolutionary step by which a lineage can explore more distant regions of the fitness landscape than it can by single small-effect mutations (91). We define “large-effect” as mutations that cause specific, replicable, and major changes in flowering time and have been extensively experimentally validated (Table 1 and ). The majority of our mutations were loss-of-function, ascertained at the gene expression and protein levels (). This is realistic for alleles segregating in natural populations of A. thaliana that underlie adaptive phenotypic variation in metabolic, defensive, morphometric, and phenological traits (14–16, 92–100). We also included a single instance of gain-of-function through introgression of a strong functional allele to replace a nonfunctional variant of the gene FRIGIDA in the reference ecotype (101). We do not speculate about the precise nature of DNA sequence-level change since there are many potential paths to large effect, but rather we seek to understand how large-effect mutations interact with the environment to influence selection on phenotypes. In order to disentangle the effects of correlated environmental signals like day length and temperature, we selected field sites in a latitudinal gradient, in which day length and temperature covary; a longitudinal gradient, in which temperature varies independent of day length; and a seasonal gradient, in which plantings were repeated within a site for several seasons () (55, 58, 102, 103). We report morphometric, phenological, and fecundity data from this common garden experiment to answer 3 main questions:

Table 1.

Description of genes mutant in this experiment with their classical pathway designation

Gene	Pathway	Effect on flowering time
FVE*	Autonomous	Pos. regulator that suppresses FLC expression
HUA2*,†	Autonomous	Neg. regulator that enhances late flowering in FRI functionals
LD*	Autonomous	Pos. regulator that suppresses FLC expression
FT*	Integrator	Florigen integrating all pathways to activate floral meristem genes
TFL2*	Integrator	Floral meristem identity; maintains high FLC despite brief warmth
GAI*	Hormone	Pos. regulator by up-regulating FT and gibberellin synthesis
SPY*	Hormone	Neg. regulator downstream of GAI that suppresses gibberellin
CO*	Photoperiod	Sensor of inductive long days, downstream of PHYs
CRY1*	Photoperiod	Pos. regulator of CO; perceives blue light of long days
GI*	Photoperiod	Integrates circadian clock information to perceive long days
PHYA*	Photoperiod	Pos. regulator of CO; perceives red light of long days
PHYB*	Photoperiod	Neg. regulator of CO; perceives red:far-red ratio of light of long day
PHYD*	Photoperiod	Neg. regulator of CO; obligately interacts with PHYB
PHYE*	Photoperiod	Neg. regulator of CO; obligately interacts with PHYB
FRL†	Vernalization	Constituent of FRI activation complex; up-regulates FLC
FRI‡	Vernalization	Activates FLC, conferring vernalization requirement
FLC†	Vernalization	FT repressor that is epigenetically down-regulated by vernalization
VIN3†	Vernalization	Required for vernalization response in FRI:FLC functionals

Pathways are not exclusive since some genes act in multiple pathways simultaneously. Some mutant genotypes harbored multiple mutations. For example, phyabde is a quadruple mutant with LoF in PHYA, PHYB, PHYD, and PHYE. Neg., negative; Pos., positive.

Mutations in this gene were not combined with any natural allele.

Mutations in this gene were induced in a background with a functional FRI allele introgressed from the Sf-2 ecotype.

Natural alleles only by introgressing a functional Sf-2 version without induced mutation.

How do large-effect mutations in flowering time genes affect phenology, morphology, and fitness in different natural environments?

How do phenotype–fitness and genotype–phenotype–fitness landscapes change with the environment?

How does gene position within a regulatory network modulate mutational effects on the genotype–phenotype–fitness map?

Results

Here, we analyze subsets of the genotype–phenotype–fitness map across environments and then combine these components. First, we examine the genotype–phenotype landscape by quantifying genotype-phenology and genotype-branching relationships in each environment. Second, we examine the phenotype–fitness landscape by quantifying phenology-fitness and branching-fitness relationships in each environment. Finally, we synthesize these across environments to describe the environmental dependence of mutational effects on phenology, branching, and fitness. Large-effect mutations in environmental signaling pathways are especially useful in these environmental comparisons because plants harboring these mutations are unable to respond to specific environmental inputs. By comparing their phenotypes to those of the same background genotype without the mutation, these mutants can be used as phytometers to gauge the relative importance of specific environmental cues for determining phenology and branching (103).

Phenology.

To describe how these mutations affected phenological traits relative to ecotype background across 8 natural environments, we used a 2-part hierarchical clustering algorithm to uncover environmental and mutant trends (Fig. 1). A linear model detected significant heterogeneity in mutant bolting times measured in photothermal units (BPTUs) among plantings, among genotypes, and among multifactorial combinations of planting and genotype (). In particular, we identified genotype × environment interactions, with the majority of mutants showing magnitude differences in mutant effects in different environments and several showing both magnitude and directional differences (). Together, these results quantify the environmental contingency of genetic effects in determining bolting time among these different environments.

Fig. 1.

Heat map of least-square means of accumulated photothermal units to bolting (BPTUs) in mutants relative to ecotype background on a log2 scale, centered within plantings, visualizing two-step hierarchical cluster by a Euclidian distance, average-based algorithm for both genotype (rows) and planting (columns). All mutant pathways are represented, although not all genotypes, since some were not planted in all 8 sites and seasons. The first word of the row identifiers shows which pathway was manipulated in a mutant, as defined by FLOR-ID (51). Lowercase gene names indicate diminished function alleles, and uppercase, functional. Colons between pathways or genes indicate multiple genetic manipulations within a line, not gene fusions. Col and Ler indicate each line’s ecotype background, Col-0 and Ler-1, respectively. Genotypes with an induced mutation combined with a functional FRIGIDA (denoted by “FRI”) were relativized against FRI Col instead of Col-0.

Specifically, we identified 3 distinct mutant classes: late-bolting, early-bolting, and sign plasticity clusters (Fig. 1). The late-bolting cluster was dominated by mutations in the vernalization and autonomous pathways known to enhance the expression of the floral repressor FLOWERING LOCUS C (FLC) or to impair its down-regulation (i.e., FRIGIDA:vernalization insensitive 3 [FRI:vin3], FRI, fve, ld), as well as mutants impaired in photoperiod response and signal integration (i.e., ft, cry1). The early-bolting cluster showed a mix of pathways. Many upstream photoperiod mutants accelerated bolting, the most consistent of which were deficient in PHYTOCHROME B (PHYB) in the Ler background. The hormone pathway mutant spy showed accelerated phenology, and gibberellic acid insensitive (gai) did not. Acceleration of bolting in the spy mutant was expected because functional SPY is a repressor of inductive gibberellins (104). The hormone pathway’s importance for accelerated flowering has often been emphasized in short days (105), but spy also accelerated flowering in summer and spring plantings that mainly experienced long days. The vernalization pathway mutations fri, hua, and flc accelerated bolting in the Col-FRI background. These loss-of-function mutations all reduce expression of FLC, allowing flowering without vernalization. Finally, the sign plasticity group consisted of genotypes that accelerated flowering in some environments and delayed flowering in others. Genotypes showing this sign plasticity most strongly included 4 mutants deficient in phytochromes (phyb, phyabe, phyad, phyabe) which are involved in photoperiod and ambient temperature sensing (106–108) as well as the photoperiod integrators gi and co (). A detailed discussion of gene-specific results is presented in the .

Branching.

Flowering time mutants pleiotropically expanded variation in branching in both ecotypic backgrounds, and, similarly to bolting time, mutant branching varied among seasons, sites, and genotypes (). Most late-bolting mutants, especially in the vernalization pathway, produced more branches, and most early-bolting mutants in the photoperiod and hormone pathways produced fewer branches than their parental ecotypes. Specifically, gai had fewer branches than its ecotype background, and the spy mutant consistently had the fewest branches across plantings (). Several phytochrome mutants had many fewer branches than either the parental ecotypes or most other mutants, consistent with the constitutive shade avoidance phenotype expected with loss of photoreceptor function. Combined with their bolting time effects, it is also possible that these mutants had fewer branches because they flowered rapidly and therefore developed fewer nodes from which branches could emerge. Several mutants also exhibited sign plasticity for branch production, including the phytochrome mutants phya, phyd, phyad, and phybe, and the vernalization mutant FRI:vin3. Flowering time mutants generated in the Col-0 background had more branches than those in Ler-1 in most plantings (). This is likely because Ler-1 harbors a loss-of-function allele of ERECTA, a positive regulator of aerial organogenesis (109), and has fewer branches than Col-0.

Univariate Selection.

The fitness effects of flowering time mutations differed in sign and magnitude across the 5 plantings for which fecundity data were collected: Norwich fall, Norwich spring, Norwich summer, Halle fall, and Valencia fall (). In particular, nearly half of all mutant genotypes switched between positive or negative relative fitness among environments (). We quantified phenotype–fitness landscapes separately for the phenology and morphometric traits for each planting. Univariate directional selection for earlier bolting and flowering occurred in 4 plantings while delayed flowering was favored in the Norwich summer planting (Table 2). Furthermore, selection for phenology measured in calendar time was stronger than for photothermal time (), which we used to scale environmentally dependent phenologies across sites (). We did not detect significant stabilizing selection in most of the plantings. Significant negative quadratic selection was detected for calendar days to bolting in Valencia, which indicated monotonic selection toward a theoretical maximum on the selection surface for an extremely early bolting time that did not fall within the range expressed by any genotype (Table 2). Directional selection consistently favored increased branching in all sites, with no evidence of stabilizing selection on branching traits (Table 2 and ).

Table 2.

Univariate selection coefficients for A. thaliana traits in 5 field environments, analogous to partial derivatives of polynomial regression techniques

Planting	DTB				BPTU				Branching
	β	SE	γ	SE	β	SE	γ	SE	β	SE	γ	SE
Halle fall	−0.144	0.023	−0.517*	0.062	−0.083	0.07	−4.222	3.180	0.766*	0.067	0.174	0.380
Norwich fall	0.204	0.13	0.021	0.08	0.192	0.092	10.014	11.045	1.136*	0.157	0.620	0.846
Norwich spring	−0.350*	0.087	−0.001	0.043	−0.270*	0.050	0.000	8.310	0.789*	0.044	−0.106	0.165
Norwich summer	0.151*	0.040	−0.071	0.038	0.159*	0.042	−6.998	6.013	0.530*	0.041	0.085	0.169
Valencia fall	−0.214*	0.060	−0.230*	0.029	−0.180*	0.057	−17.385	15.921	0.740*	0.082	0.044	0.323

β, analogous to the directional selection coefficient; BPTU, accumulated photothermal units to bolting; Branching, total branch number; DTB, calendar days to bolting; γ, analogous to the stabilizing or disruptive coefficient; SE, numerically approximated SE for the term immediately to the left.

Bolded estimates with asterisks represent significance from Bonferroni-corrected P values <0.01.

Multivariate Selection.

Taken together, univariate selection gradients indicated that fast-developing and highly branched phenotypes were generally favored. One possible genetic mechanism for reaching this optimal fitness zone is gibberellin’s simultaneous promotion of a more branched architecture and fast development; however, high gibberellin activity also suppresses flower formation early in development, which could decrease fecundity. To determine how bivariate selection on phenology and branching occurred, we modeled both in a multivariate generalized additive model (GAM) (Table 3). This showed that directional selection for greater branching was strong in all sites, but direct selection for phenology was significant in only 2 sites: Later flowering was favored in Norwich summer, and earlier flowering in Halle fall. No γ terms, which denote stabilizing, disruptive, or correlative selection, were significant. Together, this implicates the role of indirect selection on phenology while selection operated more strongly on branching (110). However, consistent with previous findings of high heritability for phenology traits and low heritability for fitness estimates (111), the ratio of genetic variation to phenotypic variation of morphometric and fitness traits (0.27 to 0.5) was much lower than phenological traits (>0.7) (). Since these mutants were selected for their effects on flowering time, phenology’s stronger genetic basis was expected but allows us to gauge the relative contribution of these underlying mutations in determining phenotype and fitness values.

Table 3.

Bivariate selection gradient analysis

Term	Halle fall		Norwich fall		Norwich spring		Norwich summer		Valencia fall
	Est.	SE	Est.	SE	Est.	SE	Est.	SE	Est.	SE
βBPTU	−0.13*	0.10	0.33	0.17	−0.01	0.08	0.20*	0.07	−0.02	0.09
βbranch	0.72*	0.05	1.15*	0.11	0.77*	0.09	0.53*	0.04	0.71*	0.08
γBPTU	−9.10	120	5.84	146	−30.2	90.8	−60.1	100	82.2	114
γbranch	0.06	0.31	0.90	1.11	−0.00	0.18	0.06	0.16	0.18	0.22
γBPTU,branch	−0.11	0.14	0.28	0.31	0.00	0.09	0.10	0.06	−0.07	0.19

β, analogous to the directional selection coefficient; BPTU, accumulated photothermal units to bolting; branch, total number of branches; Est., estimate of the coefficient; γ, analogous to the stabilizing or disruptive coefficient; SE, numerically approximated SE.

Bolded estimates with asterisks represent significance from Bonferroni-corrected P values <0.05.

Genotype–Phenotype–Fitness Landscapes.

Modeling genotype means in 2D fitness landscapes revealed remarkable variation in the genotype–fitness map among environments (Fig. 2). Ecotype background fitness relative to mutant fitness changed dramatically among environments. Although generally more fit than mutants, our reference ecotypes (Col-0 and Ler-1) showed lower relative fitness than specific mutants in certain environments. For example, in Valencia and Halle, the Col-0 ecotype was more fit than all of its mutants except one, but, in Norwich fall, both ecotypes were less fit than the majority of their mutants. Together, these patterns show that mutants achieved higher fitness through alternate phenotypic paths. For example, some mutants were invariant in phenology in certain environments but were more highly branched whereas others bolted in the direction of higher fitness but were invariant in branching. The mutations underlying these phenotypic shifts were adaptive in conditions that shifted phenotypic expression toward fitness peaks.

Fig. 2.

Fitness landscapes from generalized additive models for accumulated photothermal units bolting and total branch number. Contour line labels show fitness in seed proxy units. BPTU refers to accumulated photothermal units to bolting. Points show line averages where lines are genotypes bulked under the same maternal conditions. “Col” refers to the Col-0 ecotype average; “Ler” to the Ler-1 ecotype average; “ColFRI” refers to the genotype with a functional version of FRI introgressed from the Sf-2 ecotype. Grey vector lines represent mutations induced in the Col and Ler ecotypes, and blue vector lines represent mutations induced in the FRI (Col) genotype.

Network Effect of Mutational Perturbation.

We were interested in how gene network connectivity affected mutational perturbations on phenotypes independent of the pathway in which a mutation was induced. To answer this question, we first assessed how connected the genes in a mutant line were, summing these connections for multiple mutants. We then regressed this metric of connectivity against relativized phenotypic perturbation. Although only 3 of 15 regressions were significant after Bonferroni correction, all showed negative relationships between connectivity and phenotypic perturbation across all 3 traits and all 5 environments (Fig. 3). Overall, the magnitude of bolting time change was considerably dampened relative to branching or fitness. This indicates that bolting time is relatively robust to genetic perturbation in flowering time networks, but that branching and fitness are much less robust. In particular, mutations in more highly connected genes led to consistent and, in at least one environment, significant decreases in fitness.

Fig. 3.

Relationship between network connectedness and mutant phenotypic shift relative to ecotype background. BPTU is for accumulated photothermal units to bolting; branching is the total number of branches; and fitness is seed proxy number. Shaded areas show 95% confidence intervals; panels with asterisks denote significant linear regressions (P < 0.05) after Bonferroni correction.

Discussion

A major goal in evolutionary biology is to understand the relationship among mutations, traits, and fitness in natural environments. Directly inferring these relationships is frequently constrained by limited knowledge of genetic mechanisms in ecological systems and by the difficulty of testing mutants in multiple natural environments. However, by leveraging the well-characterized genetic architecture of flowering time and branching in the model organism A. thaliana, we investigated how mutations in multiple flowering time pathways altered the distribution of A. thaliana genotypes across fitness landscapes, the traits that composed these landscapes, and how they depended upon the environment a plant experienced.

Phenology, branching, and fitness all expressed both magnitude and sign plasticity to growth environment (Fig. 2 and ). The largest number of lines showed clear sign plasticity for fitness and the fewest for phenology, evidence that fitness was more environmentally labile than phenology. Overall, mutants were more frequently less fit than their parental backgrounds. However, for several mutants, we identified specific natural environments in which they actually had greater fitness. There are 2 explanations for this finding of greater mutant fitness only in certain environments: 1) Mutations’ phenotypic effects might differ between environments, and 2) the underlying selection surface might change. Here, we detected evidence for both mechanisms since certain environments selected for different directions and magnitudes of trait values (evidence of changing selection surfaces) and mutants’ phenological and branching effects differed among plantings (evidence of different mutant effects). Thus, we explain our patterns of conditional fitness advantage in terms of the dramatically different positions in phenotypic and fitness space that the mutants occupied upon exposure to new environments, as well as of those environments’ different selection regimes. This finding points to the importance of empirical tests of fitness in the field. While phenotypic effects may be predictable under limited laboratory conditions, the ramification of genetic perturbation on fitness in nature, mediated by environmentally influenced phenotypes, is highly environmentally contingent.

Mutation Effects on Traits.

Natural environments are far more complex than the controlled laboratory conditions in which most mutations are assessed so we might expect the unexpected when observing mutants in the wild. For example, for A. thaliana flowering time, recent empirical work showed that realistic daily amplitudes in temperature considerably decreased the delay caused by introgressing a strong FRI allele to a nonfunctional background compared to the extreme delay that occurred in constant laboratory temperatures (112). Many of our experimental mutants altered flowering in the direction we expected, but others showed complex and unexpected patterns and were likely responding to multivariate cues in their natural environments that have not yet been tested for these genes. For example, the autonomous pathway mutants lumindependens (ld) and fve were constitutively late-bolting (Figs. 1 and 2). LD and FVE down-regulate the floral repressor FLC in an environmentally insensitive manner (113–115) so we expected mutations in those genes to delay flowering in all sites. However, HUA2 is also an activator of FLC; the hua2 mutant in the low FLC ecotype background Col caused acceleration but, in the high FLC Col FRI background, switched between slight delays in Cologne, Halle, and Valencia falls and acceleration in the other plantings. We speculate that this unexpected sign plasticity may result from either 1) the ecotype’s unmutagenized background allele being responsive to an uncharacterized or multiple environmental inputs or 2) the ecotype background allele epistatically regulates an environmentally responsive module of the flowering time network. The former possibility is especially relevant for complex environments since flowering time genes have been shown to respond to a diverse set of environmental signals. For example, CO has classically been considered a photoperiod pathway gene whose primary role for flowering time is the perception of long days (114, 116, 117). However, recent work has shown that, in addition to this role, CO also integrates ambient temperature information to promote flowering (107, 118–120). The co mutant in the Ler background was delayed in the long-day plantings of summer and spring, which was likely influenced by a combination of the diminished perception of inductive photoperiods and high ambient temperatures. Remarkably, co in the Col background reversed this trend by accelerating bolting in most plantings, a dramatic demonstration of epistasis with genetic background. Thus, the effect of a mutation might change among environments as the relative importance of its roles as a photoperiod- and temperature-sensor also changes.

Fitness Landscapes.

Environmental variation in the phenotypic effects of mutations does not completely explain why the relative fitness among mutants varies across environments. Variation in the strength and direction of multivariate selection was also critically important. This is exemplified by the autonomous mutant ld, which was consistently late-bolting but switched between relative fitness increases and decreases among environments (). Thus, the topography of the selection surface shifted so that, although broad phenotypic trends remained constant, relative fitness among mutants changed. Fitness is the most composite phenotype expressed by an organism, integrating all preceding traits and developmental contingencies (111, 121). Thus, it is expected that growing conditions could perturb critical junctures during a plant’s developmental program, resulting in different fitness outcomes.

Although these mutants were chosen to expand phenological variation, they also pleiotropically expanded morphological variation in branching. Branching proved to be under stronger selection than phenology (Table 3) although both appeared to be under selection when analyzed in a univariate framework (Table 2). Predictions of responses to selection often center on highly heritable traits like phenology; however, in this experiment, phenology was under weaker selection than branching, which was much less heritable. Since a trait’s response to selection is a function of both its heritability and the strength of selection, phenology may show a marked response to selection due to its greater genetic basis, but this may represent an indirect evolutionary response to stronger selection on branching. Together, this points to the importance of considering the pleiotropy of mutations’ phenotypic effects since they may affect suites of traits with different selective potentials.

Gene Regulatory Networks.

Flowering time genes have been conceptually organized into distinct pathways such as the photoperiod and vernalization pathways (51, 122, 123). This conceptual organization is a consequence of the cue-specific manner in which the functional effects of flowering time genes are initially characterized (124, 125). However, as the interactome has been elucidated by protein–protein and genetic interaction screens, the complexity of gene regulatory networks has grown (126, 127). Indeed, recent work has shown that genes often act outside of their canonical pathways to respond to multiple environmental variables (51). For example, several “photoperiod genes” have also been shown to act as ambient temperature thermosensors (106, 119, 128). Here, we challenged mutants with natural environments that integrated multiple environmental influences. When examining phenotype deviations due to mutations, we found evidence that, indeed, more highly connected genes caused greater changes in phenotypes and greater reductions in fitness, indicating that their connections within the regulatory network were nonredundant (Fig. 3). However, this relationship was not significant for all environments or traits. This indicates that some environments exposed a gene’s connectedness more than others, possibly due to responses to different multivariate environmental inputs.

Conditionally Adaptive Mutations in Genotype–Phenotype–Fitness Landscapes.

Our results point to genotype × environment interaction as an important explanation for recent findings of nonnegligible frequencies of large-effect mutations segregating in natural populations of many species (14–16, 92–100). Contrary to the predictions of mutation-selection balance, which emphasizes inefficient negative selection to explain why mutations linger in a population (129, 130), we found evidence for balancing selection by environmentally dependent fitness. In several cases, the same flowering time pathway mutation had opposite effects on fitness in different environments, and environmentally dependent positive or negative selection acted upon conditionally adaptive traits that ramified from these specific mutations. As these targeted flowering time mutations were grown in different environments, some formed a path to high fitness through phenology, some through branching, and some through a combination of both. Phenological deviations in mutants were more predictable than branching or fitness but still differed among environments. The branching path was under stronger selection, but the phenological path had a stronger genetic basis. Thus, adaptation to natural selection on large-effect mutations is a contingent process that must be measured across multiple settings to gauge the pleiotropic adaptive value of underlying mutations.

Methods

Plant Materials.

In 2006 and 2007, the mutants were planted as part of a larger field experiment in common gardens in 5 locations (Norwich, United Kingdom; Halle, Germany; Cologne, Germany; Valencia, Spain; and Oulu, Finland) and several seasons (summer, spring, and fall) (55, 102, 103, 106, 131). Mutant genes, their primary flowering time pathway designation, and relevant environmental or genetic interactions are summarized in Table 1. All mutants were nontransgenic, induced by radiation or chemical mutagenesis, identified by screens for major effects on flowering tie and experimentally confirmed (). The single exception was the introgression of a functional FRI allele from the Sf-2 ecotype to complement loss-of-function alleles of this gene in Col-0; this allowed tests of the effect of induced mutations in backgrounds with a functional vernalization pathway (101). Seeds were bulked and stratified under common conditions (103). Germination occurred under natural light conditions in unheated greenhouses at each field site, and seedlings were transplanted to the field synchronously with germination flushes observed in natural A. thaliana populations near each field site. Between 15 and 20 replicates for each genotype were transplanted in a randomized block design. See for further details.

Phenotypic and Fitness Measurements.

We report 12 morphometric, phenological, and fitness traits (). Together with branching’s partial independence as a developmental module and its mechanistic interactions with flowering time and fitness (see Introduction), we focused on branching as a morphological trait; bolting time as a phenological trait; and fecundity as a fitness metric. Days to bolting was the number of days from transplant into the field until an inflorescence shoot was visible at the center of a rosette. To convert these measurements to photothermal time, we scaled by temperature greater than 3 °C during daylight, resulting in photothermal units to bolting (BPTUs) as in Wilczek et al. (103). Total branch number was the total number of primary, secondary, and cauline branches longer than 1 cm. We estimated fitness as the number of seeds produced by a plant in seed proxy units. These were calculated as a plant’s silique number multiplied by an average-length silique’s distance from its valve’s base to apex. Seed set is a direct measurement of fitness, and, since seeds are minute and numerous in each silique, silique length is a practical and reliable proxy for the number of seeds (132, 133). Fitness was scored as 0 for plants that died before producing siliques.

Statistical Methods.

All analyses were performed in R version 3.0.2 (134). To estimate the ratio of genetic variation to phenotypic variation for a trait, we used ANOVA procedures to calculate broad-sense heritability. In order to identify phenotype classes responding to genetic perturbation and to growth environment, we hierarchically clustered phenotypes by genotype and planting, using average, Euclidian-based clustering. To scale phenotypes for visualization in heat maps, we relativized mutant least-square means against within-planting least-square means of 1) Col-0 if the mutant was generated in a Col-0 background; 2) Ler-1 if in a Ler-1 background; or 3) Col-FRI if the mutation was in the Col-0 background with an introgressed functional Sf-2 FRI allele. Relativized traits were log2 transformed for decrease–increase symmetry, centering within plantings. The identified clusters were robust to clustering method. For phenological traits, calendar and photothermal time clusters were similar ().

To construct fitness landscapes for the 5 plantings for which fitness and growth metrics were taken, we first fit a generalized additive model (GAM) under a Gaussian distribution for each planting. To test for correlated selection on these traits, we included the interaction term BPTU × total branches. Since we wished to simultaneously estimate the magnitude, direction, and curvature of selection gradients, as well as to perform hypothesis testing on whether this selection was significant, we assessed the significance of GAM coefficient estimates with case bootstrapping procedures (32, 33, 135). This analysis combined the biological interpretability of traditional regression techniques with splines’ ability to explain complex surfaces, producing GAM coefficients that are analogous to partial derivatives of the directional selection coefficient (β) and stabilizing selection coefficient (γ) (136). P values were adjusted for multiple testing by Bonferroni corrections.

In order to quantify the network position of the mutant genes, we manually curated nonredundant genetic and physical interactions deposited in the gene network databases BioGRID 3.5.167 (137) and IntAct (European Molecular Biology Laboratory–European Bioinformatics Institute) (138). We counted each mutant gene’s experimentally validated physical or genetic interaction as a connection, and the sum of these connections as its connectedness. We considered these interactions to be edges in our network analysis, and the mutant genes as nodes. In this analysis, we regressed connectedness against relative phenotype deviation from ecotype background, reporting Bonferroni-corrected results.