Early experiences mediate distinct adult gene expression and reproductive programs in Caenorhabditis elegans

Environmental stress during early development in animals can have profound effects on adult phenotypes via programmed changes in gene expression. Using the nematode C. elegans, we demonstrated previously that adults retain a cellular memory of their developmental experience that is manifested by differences in gene expression and life history traits; however, the sophistication of this system in response to different environmental stresses, and how it dictates phenotypic plasticity in adults that contribute to increased fitness in response to distinct environmental challenges, was unknown. Using transcriptional profiling, we show here that C. elegans adults indeed retain distinct cellular memories of different environmental conditions. We identified approximately 500 genes in adults that entered dauer due to starvation that exhibit significant opposite (“seesaw”) transcriptional phenotypes compared to adults that entered dauer due to crowding, and are distinct from animals that bypassed dauer. Moreover, we show that two-thirds of the genes in the genome experience a 2-fold or greater seesaw trend in gene expression, and based upon the direction of change, are enriched in large, tightly linked regions on different chromosomes. Importantly, these transcriptional programs correspond to significant changes in brood size depending on the experienced stress. In addition, we demonstrate that while the observed seesaw gene expression changes occur in both somatic and germline tissue, only starvation-induced changes require a functional GLP-4 protein necessary for germline development, and both programs require the Argonaute CSR-1. Thus, our results suggest that signaling between the soma and the germ line can generate phenotypic plasticity as a result of early environmental experience, and likely contribute to increased fitness in adverse conditions and the evolution of the C. elegans genome.


Introduction
Phenotypic plasticity in response to environmental stress is a critical component of organismal fitness. Environmentally-induced phenotypic variation is thought to result, in part, from programmed changes in gene expression inherited through cell divisions or transgenerationally via epigenetic mechanisms such as DNA methylation, histone modifications, and non-coding RNAs [1]. For example, in nematodes, Drosophila, and humans, nutritional status during early development can modulate the longevity of subsequent generations in the absence of the original stimulus [2][3][4]. In C. elegans, the observed increased longevity is dependent upon the inheritance of starvation-induced non-coding small RNAs, likely imported into the germ line from the soma [3]. While this example is intriguing, an important and unresolved question is to what extent does environmentally-induced phenotypic plasticity mediated by changes in epigenetic marks result in adaptive variation of traits that is favored by natural selection [5,6].
The free-living nematode C. elegans is an excellent model system to investigate the molecular mechanisms regulating environmentally induced phenotypic plasticity because its developmental trajectory is dependent on the environmental conditions experienced early in life. If conditions are favorable, C. elegans undergo continuous development consisting of four larval stages (L1-L4) followed by reproductive adulthood [7]. When faced with environmental stress (e.g. over-population, low food supply, or elevated temperatures), L1 larvae initiate an alternative diapause stage named dauer. Dauers are developmentally arrested, non-feeding, nonaging larva that exit diapause only if environmental conditions are favorable [8]. Natural populations of C. elegans experience a "boom and bust" reproductive strategy whereby they exist primarily as dauers for stress resistance and geographic dispersal but reproduce rapidly when food is available [9]. A potential prerequisite for the evolution of mechanisms that modulate adult phenotypes in response to specific stressors is that early environmental conditions are predictive of future conditions. The sensitivity of this developmental system to different dauer-inducing conditions, and how this may contribute to distinct phenotypic trajectories, is unknown. Our previous work has shown that postdauer (PD) adults that experienced crowding early in development exhibit changes in gene expression, genome-wide chromatin states, small RNA populations, and life history traits compared to isogenic animals that experienced continuous development (controls, CON) [10][11][12]. Thus, numerous molecular cues may have the potential to propagate information regarding early-life experiences to modulate adult developmental outcomes.
Over a century ago, August Weismann proposed in the germ-plasm theory of heredity that only germ cells, and not somatic cells, could pass heredity information [13]. Increasingly, studies have challenged the "Weismann Barrier" by demonstrating that the passage of non-coding RNAs between somatic tissues and the germ line can result in transgenerational inheritance [14]. For example, a recent report showed that C. elegans double-stranded RNA (dsRNA) generated in neurons and transported to the germ line resulted in transgenerational silencing that is dependent upon the main systemic RNA interference (RNAi) effector, SID-1 (WBGene 00004795) [15]. Our previous work has shown that RNAi-pathways are required in different subsets of neurons for dauer formation in response to distinct environmental stresses, as well as the resulting reproductive plasticity observed between control and postdauer adults [11,16]. These observations raised the intriguing possibility that postdauer animals that experienced different early life stresses retain distinct molecular signatures mediated by non-coding RNA signals.
Here, we investigate the effect that early environmental history has on modulating phenotypic plasticity in adults. We show that postdauer adults exposed to starvation (Stv) early in life exhibit distinct gene expression profiles and reproductive phenotypes when compared to postdauer adults that experienced crowding or high pheromone (Phe). These differences are highlighted by a set of "seesaw" genes that are oppositely regulated between postdauer and control adults depending on the experienced environmental stress. In addition, we provide evidence that the significant seesaw pattern of gene expression is due to RNAi-dependent fluctuations of gene expression across whole chromosomes, resulting in transcriptome-wide seesaw trends of gene expression changes in response to stress. Moreover, our results indicate that the distinct gene expression profiles in starvation versus pheromone conditions are dependent on germline-or somatic-generated signals, respectively, that potentially move between tissue types and mediate changes in brood size. Our results suggest a model where crosstalk between the soma and the germ line governs the mRNA transcriptome and reproductive plasticity following distinct environmental histories. Finally, we provide evidence that the relationship between early life environmental stress and distinct postdauer phenotypes has contributed to the evolution of C. elegans genome organization.

Adult gene expression is dependent on environmental history
To test the hypothesis that postdauer adult phenotypes are dependent upon the dauer-inducing stress experienced early in development, we conducted RNA-Seq using wild-type (WT) control (CON Stv ) and postdauer adults that passed through the dauer stage as a result of having experienced starvation (PD Stv ) (S1A-S1D and S2 Figs). Genes with significant changes in mRNA levels due to starvation-induced passage through the dauer stage were identified by comparing CON Stv and PD Stv libraries and subjected to a false discovery rate (FDR) p-value correction of less than 5%. We identified 1,121 and 551 genes that exhibited significant upand downregulation, respectively, in wild-type PD Stv compared to wild-type CON Stv (WT Stv ) ( Fig 1A and S1 Table). To determine if PD Stv adults are distinct from PD Phe , we repeated our previous experiment [10] and conducted RNA-Seq on CON Phe and PD Phe populations grown in parallel to the starvation samples. We identified 441 and 560 genes that were significantly up-or downregulated, respectively, in wild-type PD Phe compared to wild-type CON Phe (WT Phe ) (Fig 1B and S1 Table).
Next, we compared the transcriptional changes between WT Stv and WT Phe to determine if different dauer-inducing stresses result in distinct transcriptional memories. A comparison of the 1,121 upregulated genes from WT Stv to the 441 upregulated genes identified in WT Phe identified only two genes (y38f1a.6 (WBGene00012608) and t05b11.3 (WBGene00020246)) in common. Likewise, when the 551 WT Stv downregulated genes were compared to the 560 WT Phe downregulated genes, there was no commonality (S3A Fig). However, we found that 249 (56%) of the 441 WT Phe upregulated genes were downregulated in WT Stv (WT Phe up:: WT Stv down). Similarly, 263 (47%) of the 560 downregulated WT Phe genes were upregulated in WT Stv (WT Phe down::WT Stv up) (Figs 1A-1C and S3B and S2 Table). This observation indicates that not only do PD Stv and PD Phe adults have distinct transcriptional profiles, but also that a subset of genes is oppositely regulated based on the experienced dauer-inducing stress. Due to their propensity for being up-or downregulated in postdauers compared to controls depending on the animals' environmental history, we refer to these 512 genes as the "seesaw" genes. Using random simulations, we determined that the observed numbers of seesaw genes in WT Phe and WT Stv gene sets are significantly greater than expected by chance (22.69-fold increase, p < 0.0001 and 9.38-fold increase, p < 0.0001, respectively).
To characterize the seesaw genes, we determined their expression patterns using several curation methods. First, we determined the enrichment of germline expressed seesaw genes by comparing their mRNA levels in wild-type animals to temperature-sensitive glp-4(bn2) (WBGene00006936) mutants [17]. At the restrictive temperature, glp-4(bn2) animals are deficient in germline stem cell proliferation, resulting in the lack of a functional germ line and sterility [18]. The analysis revealed that the expression level of WT Phe up::WT Stv down genes was substantially decreased in glp-4(bn2) compared to wild-type, indicating that these genes have germline enriched expression ( Fig 1D). In contrast, WT Phe down::WT Stv up genes exhibited comparable expression levels in wild-type and glp-4 germline-deficient worms (Fig 1E), suggesting that expression of these genes is enriched in somatic tissue. These trends in expression were further supported by examining the overlap of seesaw genes with tissue-enriched gene lists [19], which revealed that 73% of WT Phe up::WT Stv down seesaw genes were germlineenriched (S3C and S3D Fig and S3 Table). Lastly, we used the gene function descriptions in WormBase (WS253; [20]) consisting of published observations, GO annotations, tissue expression data, modENCODE [21], and Ensembl [22] information to further curate seesaw gene predicted function and spatial expression. Using these combined curation sources, we found an over-representation of genes associated with reproduction and embryonic development amongst the WT Phe up::WT Stv down seesaw genes [23,24] (Fig 1F and S4 Table). In contrast, the expression of WT Phe down::WT Stv up seesaw genes are enriched in the intestine and nervous system, and have putative functions associated with innate immune response or antimicrobial defense [23][24][25][26] (Fig 1G and S5 Table). Together, these results indicate that C. elegans animals maintain a cellular memory of their early life experience through the expression of sets of genes that are sensitive to environmental history and are distinct in their expression profiles and functional composition. Volcano plots representing WT Stv (PD Stv /CON Stv ) and WT Phe (PD Phe /CON Phe ) RNA-Seq data. Dotted line indicates FDR cutoff of p = 0.05; dots above the lines are significantly different in PD compared to CON for a given condition. Red and blue dots represent the WT Phe down::WT Stv up and WT Phe up::WT Stv down seesaw genes, respectively. (C) Heat map depicting the expression levels of WT Phe and WT Stv seesaw genes. The complete list of seesaw genes can be found in S2 Table. A color key and histogram for the counts and values is shown in the right panel. (D, E) mRNA levels of (D) WT Phe up::WT Stv down and (E) WT Phe down::WT Stv up seesaw genes in wild-type and glp-4(bn2) strain [17]. The x-axis represents the fragments per kilobase of transcript per million mapped reads (FPKM) and the y-axis indicates the number of genes (count) per bin. (F, G) Putative function and spatial expression of the (F) WT Phe up::WT Stv down and (G) WT Phe down::WT Stv up seesaw genes. Ã p < 0.05; ÃÃ p < 0.01; ÃÃÃ p < 0.001; ÃÃÃÃ p < 0.0001; n.s. not significant (Fisher's exact test) indicates enrichment compared to neutral expectations [20,[23][24][25][26][83][84][85]. Number of genes for each category is in parentheses. See also S1, S2, S4 and S5 Tables. https://doi.org/10.1371/journal.pgen.1007219.g001 Environmental experience and gene expression

CSR-1 is required for gene expression changes between control and postdauer adults
CSR-1 (WBGene00017641) is a C. elegans Argonaute protein that protects germline-expressed "self" transcripts from RNAi silencing through the organization of active chromatin domains and promotion of sense-oriented RNA polymerase II transcription genome-wide [27][28][29][30][31][32]. We previously showed that the CSR-1 RNAi pathway is required in early larval stages for dauer formation in response to starvation and high pheromone conditions and contributes to stable PD/CON changes in the chromatin state and gene expression for a subset of genes in adults [11,16]. In addition, we found 95% of the WT Phe up::WT Stv down seesaw genes overlapped with a previously identified list of genes targeted by CSR-1 in the germ line [29], suggesting the possibility that CSR-1 may play a prominent role in the regulation of postdauer transcriptional memory as a consequence of environmental history. In contrast, only 1.1% of genes in the WT Phe down::WT Stv up dataset have been identified as CSR-1 targets. We therefore examined if the loss of CSR-1 affected the transcriptional changes observed in WT Stv and WT Phe . Since csr-1 null mutants are sterile, we used a csr-1 hypomorph where sterility is partially rescued with a germline specific transgene [29], and performed RNA-Seq on csr-1 control adults and starvation-and pheromone-induced csr-1 postdauer adults (S1E-S1H and S2 Figs). Comparison of PD/CON gene expression levels for starvation and pheromone conditions (csr-1 Stv and csr-1 Phe ) revealed that 48 genes were significantly upregulated and 224 genes were downregulated in csr-1 Stv, while only 3 genes were significantly upregulated and 87 genes were downregulated in csr-1 Phe (Fig 2 and S1 Table). In addition, only 7 genes (abu-14 (WBGene 00004174), c47f8.7 (WBGene00008163), cut-2 (WBGene00009983), f53a9.8 (WBGene0001 8731), r02f11.1 (WBGene00019839), tts-1 (WBGene00006650), and tts-2 (WBGene00006651)) exhibited csr-1 Phe down::csr-1 Stv up seesaw patterns of gene expression, amongst which only one (f53a9.8 (WBGene00018731)) was also found to seesaw in the WT Phe down::WT Stv up dataset. None of the genes exhibited a csr-1 Phe up::csr-1 Stv down seesaw expression pattern (Figs 2 and S3E, S3F and S2 Table). These results indicate that a functional CSR-1 RNAi pathway is required for the transcriptional memory of developmental history in starvation and high pheromone dauer-inducing conditions, including for the majority of WT Phe down:: WT Stv up genes that have not been previously identified as CSR-1 targets.
We sought to verify the RNA-Seq results using qRT-PCR on a subset of genes from biologically independent samples of wild-type and csr-1 hypomorph strains. First, we examined the mRNA levels of 12 germline-specific, CSR-1-targeted genes that exhibited the WT Phe up:: Environmental experience and gene expression WT Stv down seesaw pattern. Seven (cye-1 (WBGene00000871), f45f2.10 (WBGene00018482), isw-1 (WBGene00002169), ifg-1 (WBGene00002066), cbd-1 (WBGene00010351), daz-1 (WBGene00000935), and lin-41 (WBGene00003026)) of the 12 genes (58%) were validated in wild-type animals (Figs 3, S4C and S4E); however, none of the wild-type seesaw patterns were validated in the csr-1 hypomorph (Figs 3 and S5C). For the validated genes, the abrogation of seesaw gene expression in the csr-1 hypomorph was due to multiple effects on gene expression levels in postdauers and/or controls, including: PD/CON direction of gene expression change inverting in both conditions (29%); both conditions exhibiting similar direction of change to either WT Phe (29%) or WT Stv (29%); or one or both conditions no longer exhibiting a significant change in PD/CON mRNA levels (14%) (S4C and S5C Figs).
Next, we sought to validate the mRNA levels of 12 soma-enriched, non-CSR-1-targeted genes that exhibited the WT Phe down::WT Stv up seesaw patterns. We validated 10 (ins-19 , and y51f10.7 (WBGene000 21768)) out of the 12 genes (83%) in wild-type samples (Figs 4, S4B and S4E). Consistent with our transcriptome data, the seesaw gene expression pattern for a majority of these genes were also dependent on CSR-1, despite not being previously identified as targets of CSR-1 nor being germline-enriched. Only 2 genes, hsp-16.41 and r12e2.15, retained a significant seesaw pattern in the csr-1 hypomorph; however, in both cases, the change in expression between csr-1 Phe and csr-1 Stv is opposite to the change observed in wild-type samples (S4B and S5B Figs). For the remainder of the validated genes, they showed similar disruptions in expression in csr-1 hypomorph compared to wild-type as the germline-enriched genes: the PD/CON direction of gene Environmental experience and gene expression expression inverting in both conditions (30%); both conditions exhibited trends in gene expression similar to either WT Phe (20%) or WT Stv (10%); or that one or both conditions no longer exhibited a significant change in PD/CON mRNA levels (20%). These results further confirm that CSR-1 plays a crucial role in mediating the PD/CON gene expression changes based on environmental history, regardless of whether the gene is a known CSR-1 target or not.

SID-1 requirement for seesaw gene expression
Thus far, our implication of the CSR-1 RNAi pathway in the regulation of soma-enriched, non-CSR-1 target seesaw genes suggests the possibility that signals transported between cell and tissue types could contribute to the transcriptional memory of environmental history. In C. elegans, systemic RNAi spreads dsRNA throughout the animal and requires its main effector, the dsRNA importer, SID-1 [33,34]. To ascertain whether the transport of dsRNA is a mechanism eliciting the starvation-and pheromone-induced seesaw effect, we measured the PD/CON mRNA levels of genes in the sid-1(qt9) null mutant. Similar to the csr-1 hypomorph, both germline-and soma-enriched genes failed to exhibit a seesaw pattern in sid-1 adults (Figs 3, 4, S4 and S7). When we examined whether sid-1 was required for the seesaw gene expression for a specific dauer-inducing condition, we found that 47% of the validated genes exhibited the opposite direction of change in sid-1 compared to wild-type for the pheromone condition, compared to 20% showing a similar effect in starvation condition. Thus, SID-1 primarily contributes to seesaw gene expression profiles in somatic and germline tissues due to early life history of the pheromone condition. Environmental experience and gene expression

Germline requirement for seesaw gene expression
The observation that signals from the germ line can mediate somatic gene expression levels to affect adult lifespan is well-established in C. elegans, Drosophila, and mammals [35]. Since we have shown that the systemic RNAi effector, SID-1, is playing a role in the transcriptional memory of environmental history, we interrogated whether signals exported from the germ line are necessary for the seesaw pattern in adult somatic tissue. To examine the potential role of the germ line in modulating the gene expression changes due to environmental history, we performed qRT-PCR to measure PD/CON mRNA levels for germline-and soma-enriched seesaw genes in a strain carrying the glp-4(bn2) allele, which lacks a functional germ line at the restrictive temperature [18]. As expected in animals lacking a germ line, the pheromone-and starvation-induced seesaw expression of the 12 germline-enriched genes was abolished in glp-4 (bn2) adults grown at the restrictive temperature (Figs 3 and S6C). Similarly, all but one (f55b11.4) of the 13 soma-enriched genes also showed an elimination of the seesaw pattern in glp-4(bn2) adults (S6B Fig). These results indicate that the seesaw pattern of gene expression, including genes that have enriched expression in the soma, requires a functional germ line. Again, to examine whether the germ line is required for the regulation of gene expression in specific dauer-inducing conditions, we examined whether PD/CON ratio of mRNA levels were affected for the pheromone or starvation conditions in glp-4 animals. For both germline and soma-enriched genes, we observed that 71% of the validated genes exhibited the opposite direction of change in expression in glp-4 compared to wild-type for the starvation condition, while only 24% exhibited this effect for the pheromone condition. Thus, a functional germ line is paramount for the programmed change in PD/CON mRNA levels in both the soma and germ line as a result of early life starvation. Furthermore, since SID-1 is not required for starvation-induced expression changes, these results indicate that the germline-dependent signal regulating somatic gene expression is not dsRNA.

Genome-wide seesaw patterns are chromosome dependent
The significant excess of genes exhibiting seesaw patterns of differential expression and the ability of CSR-1 to modulate expression of non-target genes led us to investigate whether the inverse expression response to distinct dauer-inducing stresses may be a genome-wide phenomenon. This analysis revealed that a large proportion of genes whose PD/CON mRNA levels were not significantly seesawing by our original, stricter criteria still exhibited opposite changes in PD/CON mRNA levels with respect to the dauer-inducing stress. A significant inverse correlation in gene expression change was observed for 67.1% of genes (12,454 out of 18,570 genes sampled in both experiments) in response to two dauer triggers (R 2 = 0.167, p < 0.0001) (Fig 5A, Q1 and Q3). Consistent with our previous analysis, most CSR-1 targets (74.3%) exhibited the WT Phe up::WT Stv down pattern of expression (Fig 5A, Q1). Moreover, when we compared the transcriptome of csr-1 Phe to csr-1 Stv , we observed a pattern distinct from wild-type, particularly for genes that upregulated in WT Phe dataset (Fig 5A and 5B, Q1 and Q4). This analysis suggests that a majority of the genes in the genome, in addition to our identified set of seesaw genes, are subject to trends of CSR-1 dependent differential regulation in response to environmental history.
A recent report found that genes that were similarly downregulated via RNAi-mediated chromatin remodeling in S. pombe during quiescence are located in clusters throughout the genome [36]. In C. elegans, thousands of protein-coding, germline-expressed genes are physically clustered in euchromatic domains that are established and maintained by the CSR-1 pathway [29,31,32]. To assess whether the physical location of C. elegans genes correlates with the observed transcriptome-wide seesaw patterns of gene regulation, we examined the location of CSR-1 target genes throughout the genome. First, examination of the chromosomal distribution of CSR-1 targets identified a highly significant enrichment on chromosomes I and III and paucity on chromosome V and sex chromosome X (Fig 5C). Second, a striking bias that parallels the distribution of CSR-1 targets was observed for the genomic distribution of genes based on their response to the starvation condition. Genes that were downregulated in WT Stv (Q1 and Q2) were overrepresented on chromosomes I and III, and genes that are upregulated in WT Stv (Q3 and Q4) are overrepresented on chromosomes V and X (Fig 5A and 5C). Expression patterns of genes in the pheromone condition did not correlate with their distribution across chromosomes beyond the seesaw relationship in expression responses observed in Q1 and Q3 (Fig 5A). Thus, our analysis shows that genes with similar trends in expression patterns in starvation condition are non-randomly distributed across chromosomes, with CSR-1 targets overrepresented on the same chromosomes as genes that are downregulated in WT Stv .
In light of the marked chromosomal bias in expression patterns of genes in response to environmental history, we next investigated whether genes with a specific expression pattern exhibited spatial organization within chromosomes. Low recombination in the center of C. elegans chromosomes has resulted in extensive linkage disequilibrium and the operation of selection at the level of large haplotype blocks [37], possibly facilitating the distribution of similarly regulated genes in clusters. A sliding window approach was used to assess the distribution of inversely regulated genes (Q1 and Q3) exhibiting at least 2-fold change in PD/CON mRNA levels for pheromone and starvation conditions. This analysis confirmed a general enrichment of Q1 genes (red lines) in central regions of chromosomes I and III, and Q3 genes (blue lines) in the central regions of chromosomes II, V, and X ( Fig 5D). Further, to define the location of CSR-1 target genes relative to each other, we employed a genome-wide clustering algorithm using 4,191 CSR-1 targets identified in the germ line [29] in order to delineate CSR-1 "clusters." We found that 73% of CSR-1 targets mapped to 507 clusters ranging from 3 to 77 genes, with most (78%) being between 3 to 10 genes long (S6 Table). The number of clusters (1.33x; p < 0.0001) and the number of clustered CSR-1 targets (1.73x; p < 0.0001) both exceeded neutral expectations as defined by randomized gene order simulations. As would be expected, there was a highly significant enrichment of CSR-1 clusters on chromosomes I and III (p = 0.0002). Since germline-expressed genes and CSR-1 targets are enriched in operons [38,39], we also examined the expression patterns of a defined set of C. elegans operons [40]. Interestingly, while 359 of the 901 known operons exhibit a WT Phe up::WT Stv down directional change in gene expression, only 91 of the 512 (18%) significant seesaw genes reside in operons, suggesting that operons alone do not account for our overarching genomic trends in gene expression. Thus, the co-localization of genes exhibiting similar seesaw trends in expression to particular chromosomes indicates that many of these loci would be simultaneously captured via genetic hitchhiking in the repeated selective sweeps that have shaped the C. elegans genome.

Early life history mediates reproductive plasticity
We next questioned whether the significant seesaw changes in mRNA levels due to environmental history result in phenotypic consequences in adult animals. Due to the overlap between germline-specific genes and WT Phe up::WT Stv down seesaw genes (Figs 1D, 1F and S3D Environmental experience and gene expression hypothesized that a reproductive phenotype could be an outcome of an animal's environmental history. To determine whether seesaw gene expression affected the number of progeny produced by wild-type postdauer hermaphrodites, we quantified the brood size of control and postdauer adults that experienced either pheromone or starvation. Consistent with our previous reports, wild-type PD Phe had an increased brood size compared to wild-type CON Phe ( Fig  6A; [10,11]). In contrast, wild-type PD Stv had a reduced brood size compared to wild-type CON Stv (Fig 6B). This indicates that the gene expression changes resulting from environmental history have significant consequences with respect to the fitness of C. elegans animals.
Since we observed that seesaw changes in gene expression are dependent on mechanisms involving CSR-1, GLP-4, and SID-1 functions (Figs 3, 4, S4, S5, S6 and S7), we asked whether the same mechanisms affected the fecundity differences observed in pheromone-or starvation-induced wild-type postdauer and control adults. Even with a reduced brood size, the number of surviving progeny produced by the csr-1 hypomorph remained slightly higher in PD Phe compared to CON Phe (Fig 6A). However, the decrease in brood size between postdauer and control adults was abrogated in csr-1 Stv (Fig 6B). We also observed that glp-4(bn2) adults grown at the permissive temperature no longer exhibited a decrease in postdauer brood size in the starvation condition (Fig 6B), but continued to show an increased number of progeny in the pheromone condition ( Fig 6A). Recently, GLP-4 was shown to be expressed in the intestine and somatic gonad in addition to the germ line [41]; thus, our observed glp-4 PD/CON brood size phenotype at the permissive temperature suggests that GLP-4 function in the intestine or somatic gonad is contributing to the regulation of the starvation program (see Discussion).
In contrast, the brood size results of the sid-1(qt9) strain were opposite to those of the glp-4 (bn2) strain. The increase in PD Phe /CON Phe brood size observed for wild-type adults was abolished in sid-1(qt9) adults ( Fig 6A); however, the decrease in PD Stv /CON Stv brood size was also observed in sid-1 adults (Fig 6B). Together, these results are in accordance with our gene expression analyses where changes in mRNA levels and the resulting phenotypic plasticity due to pheromone conditions are a result of SID-1 function in the soma, whereas starvationinduced changes are mediated by unknown signals from the germ line.
Next, we sought to further characterize the developmental differences in the germ line that could result in altered PD/CON brood sizes. Since reproduction in self-fertilizing C. elegans hermaphrodites is sperm-limited [42], we asked whether the fecundity differences in PD/CON was associated with changes in mRNA levels of genes regulating germ line mitotic proliferation and the onset of meiosis during spermatogenesis. Hermaphrodites possess two gonad arms, each of which is capped by a distal tip cell (DTC) that maintains the germline stem cell niche through GLP-1/Notch (WBGene00001609) signaling (mitotic zone, MZ) [43,44]. As cells divide in the mitotic proliferative zone, the most proximal cells begin to express the RNA binding protein, GLD-1 (WBGene00001595), which promotes entry into meiosis (transition zone, TZ) [45,46]. Since hermaphrodites produce all their sperm during the larval L4 stage, we hypothesized that modulation of these genes as a result of environmental history could potentially alter the number of sperm in hermaphrodite animals. To test our hypothesis, we first examined the expression of WT Phe up::WT Stv down seesaw gene, gld-1, using a gld-1::gfp transgene expressed in the germ line. To compare developmentally synchronized animals, we examined GFP levels of postdauer and control animals that experienced either pheromone or starvation conditions and exhibited the vulva morphology characteristic of L3, L4.1, and L4.4 larval animals [47,48], at which times the mitotic and transition zones are evident [43]. Although we were unable to validate the seesaw changes in gld-1 expression using qRT-PCR (S4C Fig), we observed that GFP levels in the germ line were significantly increased in PD Phe larva compared to CON Phe larva at all stages, consistent with our RNA-Seq results (S8A and S8B Fig and S11 Table). In contrast, we detected no significant change in GFP levels in PD Stv compared to CON Stv in animals exhibiting L3 vulval morphology, but a surprising increase in GFP levels in PD Stv compared to CON Stv for the L4.1 and L4.1 stages similar to the pheromone condition (S8A and S8C Fig and S11 Table). In Notch signaling mutants, GLD-1 levels remain stable in the transition zone [49], and high levels of GLD-1 are sufficient to drive germline stem cells into meiosis, even in the presence of Notch signaling [50]. Furthermore, we observed that the area of the gonad arms was also significantly different due to environmental history, such that postdauer gonad arm area was increased or decreased compared to controls in pheromone and starvation conditions, respectively, for all larval stages (S8B and S8C Fig and S11 Table). Together, these results suggest that the changes in GLD-1 levels we observed in Phe larva likely reflects alterations in the numbers of cells entering meiosis and not changes in GLD-1 expression in individual cells.
To further test our hypothesis, we DAPI-stained larva and counted the number of cell rows per gonad arm in control and postdauer animals that exhibited the characteristic L3, L4.1, and L4.4 vulva morphology and experienced either pheromone or starvation conditions (S9 and S10 Figs) [48,51]. If the sperm to oocyte developmental switch remains constant between postdauer and control animals [52], we would predict that postdauer larva that experienced pheromone or starvation conditions to begin germline proliferation earlier or later than control larva, respectively, resulting in altered numbers of sperm available for self-fertilization in adult hermaphrodites. Indeed, we observed a significant increase for PD Phe /CON Phe total cell rows and decrease for PD Stv /CON Stv total cell rows for all developmental stages examined (Fig 6C  and 6D and S12 Table). Interestingly, at the L3 stage, we observed a similar number of cell rows in the mitotic zone, but different numbers of cell rows in the transition zone for postdauers compared to controls in both conditions, suggesting that proliferation begins earlier or later in postdauer animals that experienced pheromone or starvation conditions, respectively (Fig 6C and 6D and S12 Table). As the animals aged and germ lines expanded, we observed that the numbers of cell rows in a particular region of the gonad, such as the transition zone or pachytene zone, varied between postdauer and control animals differently depending on the stage (L4.1 or L4.4) and environmental condition (pheromone or starvation). This result likely reflects the different mechanisms regulating pheromone and starvation gene expression changes. As an additional control, we also counted the number of cells in the spermatheca [53], which is a part of the somatic gonad, and found that the number of cells are similar for all the L4.4 populations as expected (Fig 6C and 6D and S12 Table). This result indicates that germline development, and the onset of germline proliferation, can be uncoupled from somatic gonad development, resulting in postdauer germ lines that are "older" or "younger" compared to their control counterparts. Based on our RNA-Seq data in adults (Fig 1A and  1B), these developmental trends seem to persist from the L3 stage into adulthood to result in the WT Phe up::WT Stv changes in expression of germline-enriched genes. Together, these results are consistent with the model that the onset of germline proliferation during L3 larval stages is determined by environmental and developmental history, resulting in altered sperm number and brood size in adults. Environmental experience and gene expression Furthermore, we sought to identify genes in addition to gld-1 that may contribute to altered germline development and spermatogenesis due to environmental experience. Using sperm transcriptome and proteome datasets [54], we found significant overlaps between the sperm transcriptome and proteome with the WT Phe up:: WT Stv down (p < 0.0001 and p = 0.002, respectively; two-tailed Fisher's exact test) and WT Phe down:: WT Stv up (p < 0.0001 and p < 0.0001, respectively; two-tailed Fisher's exact test) seesaw genes. GO term analyses revealed significant functional distinctions between the genes encoding protein components of sperm ("sperm genes") overlapping with the two classes of seesaw genes. Sperm genes associated with the WT Phe up:: WT Stv down seesaw genes were enriched for reproduction, genitalia development, oogenesis, and spermatogenesis, while the sperm genes overlapping with the WT Phe down:: WT Stv up seesaw genes were devoid of reproduction-related functions and were instead richly affiliated with the cuticle and collagen (S7 Table). Additional experimentation will be required to determine if the changes in sperm genes are causal to, or result from, changes in germline development in larva. In sum, these results indicate that the fecundity differences resulting from distinct life histories may have a direct relationship with the differential expression of sperm-related genes with diverging functionalities.
Taken together, we posit a model whereby different "programs" regulate global changes in gene expression leading to distinct reproductive phenotypes. In animals that experience early life high pheromone condition, the changes in PD/CON mRNA levels in the germ line and soma are dependent on SID-1 function in the soma. In contrast, in postdauer animals that experience early life starvation condition, the observed gene expression changes in the germline and the soma are not dependent on SID-1, but are instead dependent on an unidentified signal(s) from a functional germ line. Our data indicate that these two programs maintain a functional balance within the animal, such that when one program is disrupted by mutation (e.g. glp-4(bn2) mutant), the animal exhibits the gene expression and reproductive phenotype of the alternate program (Figs 3, 4, 6, S4, S5, S6 and S7). In addition, both programs are dependent on the CSR-1 RNAi pathway for these chromosomally-regulated gene expression differences to result in reproductive plasticity of adult animals (Fig 7). Moreover, we provide evidence that genes with similar trends in expression levels in response to pheromone and starvation conditions are co-localized on specific chromosomes, and may have contributed to the evolution of the C. elegans genome.

Discussion
Phenotypic plasticity in animals facilitates evolution by promoting adaptation to environmental fluctuations, colonization of novel habitats, species interaction, genetic diversification, and speciation [55,56]. C. elegans exist primarily as dauer larva in natural populations [57], suggesting that they face intense selection pressure in environments with fluctuating and often limited resources. If the early dauer-inducing environment is predictive of future conditions, we could reasonably speculate that the cellular memory of the dauer-inducing trigger could result in a gene expression program that affords postdauer adults increased fitness in a future adverse environment. For example, starvation in early life may prime a postdauer animal to curtail fecundity in order to conserve resources and minimize the consequences of another impending famine on its progeny. In contrast, over-population could reliably predict increased opportunities of outcrossing in PD Phe progeny and result in beneficial increases in the expression of reproduction-associated genes and higher brood size [58]. Thus, the phenotypic variation resulting from the pheromone and starvation programs in C. elegans may be the product of an evolutionary strategy to increase future fitness in response to past environmental challenges. Future studies will be necessary to determine whether the environmental effects on

Distinct life histories: Pheromone-and starvation-induced programs
Our results support a model where at least two distinct life history programs, the pheromone and starvation programs (Fig 7A), are orchestrated in postdauers through the exchange of signals between the germ line and the soma (Fig 7B). An intriguing question is what the candidate somatic and germline signals regulating the PD/CON programs might be. In the pheromone program, the requirement for SID-1 strongly implicates dsRNAs as mediators of gene expression changes in the soma and the germ line. Endogenous dsRNAs could potentially be processed into, or stimulate the production of, siRNAs that associate with CSR-1 in specific tissue types, analogous to how exogenous dsRNA results in target specific silencing by Argonautes such as NRDE-3 (WBGene00019862) [59]. Since one of the proposed functions of CSR-1 is the establishment and maintenance of euchromatic chromatin states associated with target genes, and we have shown that CSR-1 targets exhibiting similar gene expression patterns cluster at the chromosomal level (Fig 5D), it is reasonable to propose that CSR-1 would also modify chromatin states due to pheromone and starvation in a tissue-specific manner. Indeed, we previously showed that changes in the levels of two histone modifications associated with euchromatin in PD Phe /CON Phe are dependent on CSR-1 [11]. These altered chromatin states could potentially spread along chromosomal regions including non-CSR-1 target genes, resulting chromosomal-biased changes in gene expression genome-wide (Fig 5).
In contrast, the starvation program is independent of SID-1 but dependent on the germ line ( Fig 7B). One candidate mechanism for regulating the starvation-induced cellular memory is endocrine signaling, which plays a key role in lifespan and stress response in worms, flies, and mammals [60]. Effectors of endocrine signaling in C. elegans include the conserved FOXO transcription factor DAF-16 (WBGene00000912) and the nuclear hormone receptor DAF-12 (WBGene00000908), both of which function downstream of the insulin/IGF-1 signaling (IIS) pathway [61]. In C. elegans, Drosophila, and mice, removal or alteration of the germ line results in somatic aging phenotypes [60]. Nematodes with defective gonads exhibit extended lifespans that are dependent on DAF-16 activity in the intestine, as well as DAF-12 activity in the somatic gonad [60,62], suggesting that the regulation of a somatic phenotype is mediated by germline signaling. Recent reports provide further evidence of such crosstalk: the chemotaxis response to the odorant diacetyl in C. elegans is dependent on DAF-16 and germline proliferation [63]; cold tolerance in worms is dependent on a feedback mechanism involving IIS signaling, temperature-sensing neurons, the intestine, and sperm cells [64] while a thermosensory neuronal circuit promotes longevity at warm temperatures by engaging endocrine signaling [65]; and the transgenerational lifespan extension in a mutant strain with defective chromatin remodeling is dependent on a germ line to soma signaling mechanism modulated by DAF-12 [66]. Therefore, endocrine signaling is an attractive contender with which distinct life history trajectories program a postdauer animal in response to early life starvation.

Life history and the evolution of genome organization
Physical co-localization of genes with correlated expression patterns is widespread in eukaryotes [67] and is pronounced amongst genes expressed in the testis [68] and those encoding protein components of sperm [69]. Consistent with the efficient co-regulation of neighboring genes, we observed an excess of clustered genes in the central portions of chromosomes I, II, III, V, and X that exhibited WT Phe up::WT Stv down and WT Phe down::WT Stv up trends in gene expression changes (Fig 5A and 5D). Given the breadth of syntenic conservation between C. elegans and C. briggsae [70], it would appear likely that this genomic architecture arose in a common ancestor and may have contributed to reproductive fitness in what was likely to be an out-crossing species. Comparative genomic analyses will thus be essential in reconstructing the evolutionary history of PD reproductive investments using C. briggsae, which also shares a conserved CSR-1 RNAi pathway [38]. However, it is also noteworthy that the non-random organization of differentially expressed genes identified in this study overlaps substantively with regions of the C. elegans genome that experience limited recombination and large-scale selective sweeps [37]. Population genetic models predict that recurrent selective sweeps in linked regions would favor the establishment of adaptive "supergenes" over biologically realistic time scales [71,72]; thus, we speculate that this process may have contributed to the expansion of regions harboring genes with coherent reproductive functions during C. elegans evolution [37,73,74].
To collect PD Phe , we used an egg white plate procedure described previously [76]. To obtain PD Stv , well-fed worms were transferred to seeded NGM plates and monitored until the E. coli OP50 food was depleted and the plates were populated with dauers. See Supplemental Experimental Procedures for details.

Brood assays
Ten L4 larvae were singled onto seeded NGM plates and transferred daily onto fresh plates until egg laying ceased. Only surviving progeny were counted. At least three independent biological replicates were conducted. Statistical significance was determined using GraphPad Prism v.7.

RNA extraction, RNA-seq library preparation, and data analyses
Total RNA extraction was done using TRIzol Reagent (Life Technologies). Two biological independent RNA-Seq libraries for a strain and a condition were prepared using the NEBNext mRNA Library Prep Master Mix Set for Illumina (NEB). Data analysis was conducted on the CLC Genomics Workbench v.8.5 (Qiagen) with differential expression calculated using EdgeR [77]. Detailed procedures are described in Supplemental Experimental Procedures.

Gene ontology term analyses
GO terms were analyzed using DAVID (Database for Annotation, Visualization and Integrated Discovery) v.6.7 [78].

Monte carlo simulation
The expected frequency of seesaw genes was investigated using a randomized approach whereby sets of genes, equal in number to the observed sets of differentially expressed genes under each condition, were selected from the whole genome set, without replacement, and their overlap assessed. The expected number of overlapping "seesaw" genes and the significance of the observed numbers were directly determined from the simulated distributions, based on 10,000 simulations.

Quantitative reverse transcription PCR
qRT-PCR was done using samples collected from three independent biological replicates. S8 Table lists the primer sequences. Statistical analysis was done using GraphPad Prism v.7. Detailed procedures can be found in the Supplemental Experimental Procedures.

Genomic clustering of CSR-1 targets and direction of change
The Global Landscape Clustering (GLC) algorithm (Borziak et al., manuscript in preparation) was used to identify maximal sets co-localized genes that share a specific attribute, such as being CSR-1 targets. See Supplemental Experimental Procedures for a detailed description. Significance of clustering was directly determined from the simulated distributions, based on 10,000 simulations.

Sliding window analysis
Seesaw gene enrichment sliding window analyses were conducted on the EdgeR generated fold changes using gene position information based on the WormBase (WS235) annotations, using only coding transcripts. Sliding windows of 2.5 megabase pairs with 500 kilobase pair intervals were used. Seesaw genes were defined as those showing the inverse directional change of at least 2-fold, regardless of significance. Enrichment was calculated against the genomewide average.

Genome-wide gene expression graphs
Genome-wide gene expression graphs were generated using the log 2 transformed EdgeR fold change values. Significance of correlation between experiments was calculated based on the Pearson's correlation between samples, where degrees of freedom equals (the number of genes in the plot-1). Enrichment of directional change across chromosomes was calculated using χ 2 test with Yates' correction against the remaining chromosomes.

Germ line imaging and measurements
Detailed collection methods for BS1080 PD and CON larva in Phe and Stv conditions are described in S1 Text. For the Phe conditions, control and dauer larvae were collected using water or crude pheromone, respectively, on dauer formation plates as previously described [79; 80]. Images of larva exhibiting vulva morphology characteristic of the L3, L4.1, and L4.4 larval stages were analyzed for GLD-1::GFP expression and gonad area using ImageJ software (NIH). These same worms were next used for germ cell row counts using a standard whole worm DAPI staining protocol [81]. The size of mitotic zone, transition zone, pachytene zone was determined based on the germ cell nuclear morphology [82]. Statistical significance between CON and PD samples was determined using Student's t-test.

Accession numbers
The accession number for the high-throughput sequencing data reported in this study is GSE92954. PCA was based on the total exon reads per gene, across 21666 genes. Significance of PC loadings per sample was assessed using the FactorMi-neR dimdesc function [86]. Significance of PCs was assessed using both the Kaiser criterion and the broken stick model [87]. For both the A) wild type and B) csr-1 hypomorph PCA, only PC1 was significant using either test.   Table. qRT-PCR data for Figs 3, 4, S3, S4, S5 and S6. (XLSX) S10 Table. Raw data for brood size assays in Fig 6. (XLSX) S11 Table. Raw data for GLD-1::GFP fluorescence and gonad arm area. (XLSX) S12 Table. Raw data for germ cell row counts.