Balancing selection shapes density-dependent foraging behaviour.

The optimal foraging strategy in a given environment depends on the number of competing individuals and their behavioural strategies. Little is known about the genes and neural circuits that integrate social information into foraging decisions. Here we show that ascaroside pheromones, small glycolipids that signal population density, suppress exploratory foraging in Caenorhabditis elegans, and that heritable variation in this behaviour generates alternative foraging strategies. We find that natural C. elegans isolates differ in their sensitivity to the potent ascaroside icas#9 (IC-asc-C5). A quantitative trait locus (QTL) regulating icas#9 sensitivity includes srx-43, a G-protein-coupled icas#9 receptor that acts in the ASI class of sensory neurons to suppress exploration. Two ancient haplotypes associated with this QTL confer competitive growth advantages that depend on ascaroside secretion, its detection by srx-43 and the distribution of food. These results suggest that balancing selection at the srx-43 locus generates alternative density-dependent behaviours, fulfilling a prediction of foraging game theory.

INTRODUCTION

The benefit of a particular foraging strategy varies based on the behavior of competitors, and therefore balancing selection can favor the co-existence of multiple strategies within a species[1,2]. The pioneering example of strategic competition is natural genetic variation at the foraging (for) gene in Drosophila melanogaster larvae[3,4]. Two for alleles for active (rover) or sedentary (sitter) behavior are maintained in a population because of frequency-dependent balancing selection against larvae with the more common foraging strategy[5]. This example and others like it suggest that animals could benefit from detecting and responding to competitors in real time, and modifying foraging behavior accordingly. However, little is known about the genes and neural circuits that incorporate information about conspecifics into foraging strategies.

An opportunity to address this question is provided by the nematode C. elegans, an animal with well-characterized foraging circuits and intraspecific pheromone signaling. C. elegans foraging on bacterial food spontaneously alternates between an exploratory behavior called roaming and a less active behavior called dwelling, each of which persists for several minutes per episode[6,7]. Transitions between roaming and dwelling are regulated by distributed neuromodulatory systems that link internal cues such as nutritional status to locomotion circuits[7,8]. C. elegans senses population density using a family of secreted pheromones called ascarosides, which control the developmental decision to enter the starvation-resistant dauer larva stage[9] and regulate behaviors such as aggregation and male attraction to hermaphrodites[10,11]. Here we show that physiological levels of certain ascarosides also regulate foraging by suppressing roaming behaviors. By characterizing differences in pheromone sensitivity in natural C. elegans isolates, we identify a pheromone receptor that shapes alternative foraging strategies and affects fitness depending on the structure of the food environment.

RESULTS

Pheromones regulate foraging behavior

The effects of ascaroside pheromones on C. elegans foraging were examined by quantifying long-term exploration of a bacterial lawn by individual wild-type N2 animals (Figure 1A). To mimic the effects of high density on these isolated animals, we conducted the assay in the presence of natural pheromone extracts. The pheromones strongly suppressed exploration (Figure 1B), as did several pure synthetic ascarosides at concentrations at or below those that induced dauer larva development (Figure 1C, 1D). However, ascr#5, a potent regulator of dauer development, only weakly suppressed exploration (Figure 1D). Thus a subset of ascarosides regulates foraging behavior at biologically relevant concentrations.

Figure 1

Ascaroside pheromones suppress exploratory foraging behavior

A) Exploration assay.

B) Wild-type N2 response to crude pheromone extract, showing exploration scores and a pheromone response index, presented as mean ± SEM.

C) Structures and names of selected ascarosides.

D) N2 response to individual ascarosides, presented as mean ± SEM. ***P<0.001, **P<0.01 by ANOVA with Dunnett correction; ns, not significant.

The exploration assay is an indirect measure of the relative time C. elegans spends in roaming and dwelling states[7,8]. Quantitative behavioral analysis of video recordings showed that the potent ascarosides icas#9 or ascr#8 decreased the fraction of time spent roaming, the duration of roaming states, and locomotion speed during roaming, but did not affect dwelling states (Extended Data Figure 1).

Extended Data Figure 1

Roaming and dwelling states in the presence of ascarosides

a) Roaming and dwelling behaviors scored from video analysis. n=102-214 tracks per data point.

b and c) Cumulative distribution of roaming (b) and dwelling (c) state durations for animals in (a). ***P<0.001 by log rank test; ns, not significant.

d and e) Scatter plot of average speed and angular speed (a measure of turning rate) in 10 s intervals taken from 1.5 hour long video recordings of wild-type animals in control (d) and icas#9 (e) conditions. Roaming animals move quickly and turn infrequently compared to dwelling animals. Note bimodal distribution defining distinct behavioral states. Control = 161 tracks, icas#9 = 102 tracks.

f and g) Speed following a reversal (f) and reversal rate (g) for roaming or dwelling animals. Roaming speed is slightly slower in ascarosides (e,f). Data presented as mean ± SEM. ***P<0.001, *P<0.05 by ANOVA with Dunnett correction; ns, not significant.

A QTL for pheromone sensitivity

A variety of genetically diverse wild-type C. elegans strains responded to ascarosides with a suppression of exploration, like the control N2-like strain CX12311[12] (Figure 2A). Among the wild strains, the German strain MY14 failed to respond to 10 nM icas#9 in the exploration assay, while responding normally to ascr#2, ascr#3, and ascr#8 (Figure 2A; Extended Data Figure 2). Coupled alterations in pheromone signaling and detection can contribute to reproductive isolation during incipient speciation[13]. However, MY14 and CX12311 were found to produce similar levels of icas#9 and 16 other ascarosides (Extended Data Figure 3), indicating that the change in the icas#9 response in MY14 was independent of icas#9 production.

Figure 2

Natural genetic variation in pheromone sensitivity

A) Response of wild-type C. elegans strains to synthetic ascarosides. ***P<0.001 by ANOVA with Dunnett correction.

B) icas#9 response of 94 CX12311-MY14 recombinant inbred lines (RILs) and parental strains.

C) QTL analysis of RILs shown in B. Horizontal line denotes P<0.05 genome-wide significance threshold. LOD, log likelihood ratio.

D) Mean icas#9 response of near-isogenic lines (NILs) used to map the roam-1 locus, in MY14 and CX12311 (i) or N2 (ii,iii) backgrounds. ***P<0.001, **P<0.01 by ANOVA with Dunnett correction; ns, not significant.

All data are presented as mean ± SEM

Extended Data Figure 2

Roaming and dwelling behavior of MY14

a) Fraction of time MY14 animals spend roaming or dwelling in control, ascr#8, and icas#9 conditions. n=66-109 tracks per data point. Assays were conducted in 8% O2. Compare Extended Data Figure 1.

b and c) Cumulative distribution of roaming (b) and dwelling (c) state durations for MY14 animals scored in (a). Roaming states are significantly shorter in the presence of ascr#8 (t½ = ~150 s, vs ~220 s in controls), but they are not significantly affected by icas#9 (t ½ = ~190 s). Roaming states may also be longer at baseline in MY14 than in N2 (see Extended Data Figure 1). ***P<0.001 by log rank test; ns, not significant.

Extended Data Figure 3

Ascarosides produced by wild-type strains

LC-MS/MS analysis of ascarosides secreted by N2, CX12311, and MY14 strains grown on (a) OP50 or (b) HB101 bacteria. icas#9 is produced at similar levels by icas#9-sensitive and icas#9-resistant strains. n=2 (a) or 3 (b) culture extracts per genotype.

To determine the genetic basis of icas#9-insensitivity, 94 Recombinant Inbred Lines (RILs) were generated from intercrosses between MY14 and CX12311. A continuous distribution of icas#9 sensitivity was observed in exploration behavior of the RILs (Figure 2B) suggesting that two or more loci contribute to icas#9-sensitivity. The 94 RILs were genotyped at ~185 kb resolution across the genome by low-coverage whole-genome sequencing[14] (Supplemental Table 1). Quantitative trait locus (QTL) analysis identified a single significant QTL at genome-wide significance that accounted for 34.9% of the total variance between the RILs, roam-1 (Figure 2C). Covariate analysis failed to find additional QTLs that were either additive or interactive with roam-1 (Extended Data Figure 4).

Extended Data Figure 4

Covariate analysis of 94 RILs

Covariate analysis controlling for roam-1 genotype, testing for additive (a) or interactive (b) QTL at other loci. The horizontal line denotes the P<0.05 genome-wide significance threshold. LOD, log likelihood ratio.

The impact of roam-1 on foraging was confirmed by creating near-isogenic lines (NILs) with small genetic regions substituted between the strains. NILs in which 2.5 Mb surrounding roam-1 were reciprocally exchanged between CX12311 and MY14 were intermediate in icas#9-sensitivity compared to the parental strains (Figure 2Di). Depending on the direction of the introgression, the major-effect roam-1 locus accounted for 39-46.4% of the genetic variance between the two parental strains.

To simplify further mapping, the roam-1 region from MY14 (kyIR139) was crossed into the N2 laboratory strain. The resulting NIL (kyIR144, Figure 2Dii) facilitated further mapping that localized roam-1 to 182 kb (Figure 2Dii). High-density mapping of 2600 F2 progeny of a cross between N2 and the NIL kyIR147 yielded 12 informative recombinants in this 182 kb region that mapped the roam-1 QTL to 37 Kb on chromosome V (Figure 2Diii, Methods).

roam-1 affects the icas#9 receptor SRX-43

The 37 kb roam-1 region contains sixteen protein-coding genes, including five genes that encode predicted G protein-coupled chemoreceptors in the srx or str gene families (Figure 3A). We hypothesized that one or more of the chemoreceptors could be icas#9-receptor(s) with reduced activity in MY14, and therefore introduced N2-derived sequences overlapping the chemoreceptor genes into the roam-1MY14 NIL strain kyIR163. N2-derived transgenes covering srx-43 (bar in 3A) conferred icas#9-sensitivity to roam-1MY14 (Figure 3B), whereas transgenes with nonsense mutations disrupting the coding region of srx-43 did not (Figure 3B).

Figure 3

The roam-1 QTL includes the icas#9 receptor SRX-43

A) The roam-1 locus.

B-D) icas#9 sensitivity as affected by (B) high-copy srx-43 transgenes (C) srx-43(lf) mutations (D) single-copy integrated srx-43 transgenes.

E) srx-43::GFP translational reporter (top) and Nomarski image (bottom). Arrowhead, ASI sensory cilium. Scale bar, 10 μm.

F) ASH calcium responses in ASH::srx-43 transgenic strain.

G) icas#9 sensitivity after exchanging promoters in single-copy srx-43 transgenes.

H) Endogenous srx-43 mRNA levels.

All data are presented as mean ± SEM. ***P<0.001, **P<0.01, *P<0.05 by ANOVA with Dunnett correction or t-test; ns, not significant. Color indicates roam-1 genotype.

The function of srx-43 was examined further by characterizing loss-of-function mutations in the endogenous srx-43 gene (Methods). srx-43(lf) mutants were profoundly insensitive to icas#9 in both N2 and roam-1MY14 genetic backgrounds (Figure 3C). icas#9 insensitivity in the N2 srx-43(lf) mutant was rescued by an N2 srx-43 transgene (Figure 3C). These results indicate that srx-43 is necessary for the icas#9 response in both N2 and MY14 strains, and essential for the behavioral difference between them.

The activity of the N2 and MY14 srx-43 genes was compared by targeting a single copy of srx-43 from each strain to a defined locus using the Mos1 transposase, in an srx-43(lf) mutant so that the single-copy transgene was the sole source of srx-43. The N2 srx-43 genomic region fully rescued the icas#9 response, whereas the MY14 region did not (Figure 3D). The differential effects of single-copy transgenes indicate that MY14 srx-43 has reduced activity compared to N2 srx-43.

Reporter genes with N2 or MY14 srx-43 sequences driving GFP were expressed selectively in the ASI sensory neurons (Extended Data Figure 5A), which promote roaming behavior[8,15]. A genomic clone with GFP fused to the C-terminus of the SRX-43 protein was enriched in ASI sensory cilia, the site of sensory transduction (Figure 3E). These properties suggest that SRX-43 is a chemoreceptor. We investigated effects of icas#9 on ASI activity using in vivo calcium imaging, but failed to observe a response. This negative result is consistent with studies of dauer formation, where ascarosides regulate gene expression in ASI and not acute ASI calcium levels[16-19].

Extended Data Figure 5

Signal transduction by SRX-43

a) Expression of Psrx-43::srx-43::SL2::GFP bicistronic reporter transgenes bearing N2 (top) or MY14 (botton) srx-43 sequences. Arrows indicate cell bodies of ASI sensory neurons. Scale bar, 50 μm.

b) ASH sensory neurons are insensitive to multiple ascarosides. ASH calcium imaging with GCaMP3 in control animals that do not express srx-43 transgene, isolated as non-transgenic siblings of transgenic animals tested in Figure 3f (n=19). Ascarosides tested at 10 nM.

c) SRX-43 from MY14 confers icas#9 sensitivity on ASH neurons. Compare SRX-43 from N2 in Figure 3f.

d) icas#9 decreases daf-7::GFP expression in ASI neurons of N2 but not roam-1MY14 adults. Bars indicate mean fluorescence intensity ± SEM. *P<0.05 by ANOVA with Tukey's multiple comparisons test. n = # of animals,

e) icas#9 responses of daf-7(lf) mutants are attenuated in N2 but not in roam-1MY14 genetic backgrounds. Modified exploration assays were conducted on strains including daf-3(lf) alleles (see methods). *P<0.05, ns, not significant by t test. Data represented as mean ± SEM.

f) Time course for icas#9 response in exploration assay. Pheromone response expressed as mean ± SEM for 2, 4, 6, 10, and 14 hours following initiation of exploration assay. ***P<0.001, ns, not significant by t test with Bonferonni correction comparing squares entered in control versus 10 nm icas#9 plates. n = 12 for all time points.

To ask whether SRX-43 could be an icas#9 receptor, srx-43 cDNAs were expressed in the ASH sensory neurons, which are normally insensitive to ascarosides[17] (Extended Data Figure 5B), and ascaroside-induced calcium flux was monitored using genetically-encoded calcium indicators[17,20]. ASH neurons expressing SRX-43 responded with calcium transients to 10 nM icas#9 but not to other ascarosides or to indole (Figure 3F). Although the MY14 strain was largely insensitive to icas#9 in foraging assays, MY14 SRX-43 also detected icas#9 when expressed in ASH (Extended Data Figure 5C).

A MY14 srx-43 GFP reporter gene was expressed in ASI, but appeared weaker than the N2 GFP reporter (Extended Data Figure 5A). To ask which sequences distinguished N2 and MY14 srx-43 activity, their promoter and coding regions were exchanged and tested as Mos1-mediated Single Copy Insertion (MosSCI) srx-43 transgenes. A transgene with the N2 promoter region and MY14 coding region rescued icas#9 sensitivity in srx-43(lf) mutants but the converse did not, localizing the difference to the srx-43 promoter (Figure 3G). Quantitative measurements of endogenous srx-43 mRNA levels demonstrated that srx-43 was expressed at a five-fold lower level in roam-1MY14 than in N2 (Figure 3H). Therefore, natural variation in the activity of the srx-43 promoter between N2 and MY14 impacts srx-43 gene expression and behavioral sensitivity to icas#9.

Ascarosides promote dauer larva development in part by suppressing transcription of daf-7, which is expressed in ASI and encodes a secreted TGF-β-related peptide[19]. daf-7(lf) mutants have reduced levels of roaming, like animals treated with ascarosides[7]. Expression of a daf-7::GFP reporter was significantly reduced by icas#9 treatment, supporting a role for daf-7 as a target for icas#9, and animals bearing a daf-7 mutation were less responsive to icas#9 than controls (Extended Data Figure 5D,E). icas#9 influenced behavior only after several hours of exposure (Extended Data Figure 5F), a delay that is in agreement with slow transcriptional regulation in the daf-7 signaling pathway. Together, these results suggest that icas#9 acts as a primer pheromone[21,22] that regulates foraging via transcription and endocrine signaling.

Balancing selection at a foraging QTL

To understand the population genetics of roam-1, we examined the genomic sequence of a 20 kb region centered around srx-43 in 39 additional wild C. elegans isolates sequenced by the Million Mutation project[23]. Two discrete and highly divergent haplotypes for the roam-1 region were found, one resembling N2 and present in 34 strains, and the other resembling MY14 and present in seven strains with different geographical origins and genetic backgrounds (Figure 4A,B, Extended Data Figure 6A). The N2 and MY14 haplotypes differed at 2.64% of all positions over the 20 kb srx-43 region, 12 times the genome-wide average[24]. These data derived from Illumina sequencing that can fail to align highly divergent sequences; targeted Sanger sequenced revealed that MY14 and N2 actually differed at 19.7% of all positions in srx-43 promoter and coding regions (Extended Data Figure 6B). A phylogeny constructed for srx-43 and the most closely related genes in C. elegans, C. briggsae, and C. remanei confirmed that the divergent srx-43 alleles represent the same gene (Extended Data Figure 6C). All seven tested MY14-like strains were relatively resistant to icas#9 compared to N2-like strains (Figure 4C). These results suggest that naturally occurring icas#9 resistance is associated with a highly divergent roam-1 haplotype including srx-43.

Figure 4

Population genetics of the roam-1 locus and icas#9 sensitivity

A) A dendrogram across 41 natural isolates representing a 20 kb region surrounding srx-43. The two major roam-1 haplotypes from N2 and MY14 strains are indicated in red and blue, respectively.

B) Whole-genome dendrogram for the same strains as A, showing relationships among strains. Red and blue colors follow roam-1 haplotypes in A.

C) icas#9 responses of natural isolates with the MY14 haplotype (blue) are consistently lower than those with the N2 haplotype (red). Pheromone response shown as mean ± SEM.

Extended Data Figure 6

Alternative roam-1 allele have high sequence variability

a) The roam-1 QTL region (top). roam-1 SNPs = SNPs with respect to the N2 reference genome that are shared by JU360, MY2, MY14, ED3021, JU1171, MY16, and MY6 and not by any other strains, according to the Million Mutation Project, defining the roam-1MY14 haplotype. Other SNPs = all other SNPs with respect to the N2 reference genome found in any of the 40 wild isolates in the Million Mutation Project.

b) Polymorphisms revealed by Sanger sequencing of srx-43 promoter and coding region. Despite the high rate of polymorphism, there are only 4 nonsynonymous mutations in the MY14 coding sequence. We confirmed that the MY14 and N2 sequences are alleles of the same gene by examining sequence reads of the MY14-like strain MY23 in the CeNDR dataset (www.elegansvariation.org) and aligning each read to N2 and MY14 sequence for the srx-43 region as determined by Sanger sequencing. 7272 of the MY23(MY14) reads better matched the MY14 Sanger sequence and 4 of the reads better matched the N2 reference sequence, as expected if MY14 and N2 each bear one alternative allele of the gene.

c) Phylogeny constructed for srx-43 and related genes in C. elegans, C. briggsae, and C. remanei demonstrates that the srx-43 alleles in N2 and MY14 are closely related alleles of a single gene. Genes are color coded by species (green = C. elegans, blue = C. briggsae, orange = C. remanei). Protein sequences and gene names are from Thomas and Robertson, 2008[38].

The striking allelic divergence of the roam-1 region, combined with the observation that most genes in the interval including srx-43 have a low dN/dS ratio (Extended Data Table 1), suggests that these alternative haplotypes might be subject to balancing selection. To examine this possibility, we analyzed a database of wild strains assembled and sequenced by the Andersen lab at Northwestern University (CeNDR – www.elegansvariation.org). The MY14 haplotype was present in 21 of the 152 unique strains, and the remainder had the N2 haplotype. Both haplotypes were found globally: the MY14 haplotype was found in Europe, the United States, New Zealand, and Chile. In almost all cases in which individuals with the MY14 haplotype were isolated, individuals of the N2 haplotype were isolated from proximal environments at the same time, a distribution consistent with balancing selection.

Extended Data Table 1

dN/dS for srx-43 and other genes in the roam-1 region.

Gene	Chr.	Start	Stop	substitutions	Non-synonymous	Synonymous	dN/dS
str-231	V	16006304	16007605	36	9	27	0.097
str-233	V	16008332	16009664	85	24	61	0.118
C06C6.1	V	16010710	16012112	72	22	50	0.128
T10C6.2	V	16012537	16014226	60	24	36	0.195
str-232	V	16014509	16015752	60	20	40	0.145
srx-43	V	16017971	16019862	30	3	27	0.033
srx-44	V	16020342	16022652	25	5	20	0.071
T10C6.5	V	16023351	16024185	10	1	9	0.03
T10C6.6	V	16024326	16027454	27	2	25	0.025
T10C6.7	V	16027449	16031713	12	10	2	1.26
T10C6.16	V	16032471	16034961	28	22	6	0.943
T10C6.10	V	16036413	16038293	39	22	17	0.357
his-4	V	16040029	16040548	6	0	6	0
his-3	V	16040771	16041288	13	0	13	0
his-2	V	16041828	16042283	3	0	3	0
his-1	V	16042486	16043013	7	0	7	0
T10C6.15	V	16044296	16045874	19	9	10	0.228
ZK285.2	V	16046873	16047446	8	3	5	0.175
str-198	V	16047808	16049147	5	3	2	0.401
best-12	V	16049492	16051554	13	6	7	0.227
F14H3.3	V	16051923	16053510	32	22	10	0.522
F14H3.4	V	16053900	16055072	5	3	2	0.388
F14H3.5	V	16055806	16057054	11	10	1	2.711
F14H3.6	V	16057322	16058892	11	9	2	1.24
F14H3.15	V	16059781	16060235	11	9	2	1.37
fbxa-100	V	16060526	16062093	11	6	5	0.303
F14H3.9	V	16067949	16068763	9	7	2	0.83
cyp-35D1	V	16069238	16071301	94	27	67	0.113
nhr-176	V	16071325	16072738	15	5	10	0.127

Sequence features of the roam-1 region also matched expectations for a region under balancing selection. The region encompassing roam-1 had a relatively high Tajima's D of 1.01 that is unusual both at the genomic level (<3.4% of bins had a higher value) and in the center of chromosome V where srx-43 lies (<3.6% of bins had a higher value). A phylogenetic analysis of the 152 strains revealed that the roam-1 haplotype extends approximately 30 kb before being disrupted by recombination events (Extended Data Figure 7). Given the low outcrossing rate in wild C. elegans populations[25,26], this short roam-1 haplotype suggests co-occurrence of both haplotypes within interbreeding populations over many generations.

Extended Data Figure 7

Significant recombination between roam-1 and surrounding regions

Top, phylogenies constructed with 152 diverse wild isolates revealing differences for the region surrounding srx-43 and the regions immediately to the left and right of the 30 kb haplotype. Bottom, graph showing the number of variants and Tajima's D score calculated for 5 kb bins across a 250 kb region. The bin containing srx-43 has 250 polymorphisms and a Tajima's D of 1.01, which is high both at the genomic level (<3.4% of bins had a higher value) and for the chromosomal location of srx-43 (<3.6% of bins had a higher value).

To directly test the possibility that selection could act on roam-1, we designed competition experiments to compare the relative fitness of N2 and roam-1MY14 strains. Experiments were conducted under high-density conditions to permit the accumulation and detection of endogenously produced icas#9, and competition was applied by growing cultures past the point of starvation, i.e. with limiting food. The first competition experiments were conducted on a standard lawn of OP50 E. coli bacteria with a population founded by 20 N2 and 20 roam-1MY14 age-matched adults (Figure 5A). These conditions resulted in a growth advantage for the N2 genotype over roam-1MY14 in the first cycle of competition that continued in subsequent cycles (Figure 5B).

Figure 5

Bidirectional competitive selection at the roam-1 locus

A) Diagram of “boom-bust” competition experiments, with food depletion followed by 48 hours of starvation before transfer.

B-D) Competition on “simple lawn” showing allele ratio of DNA harvested at transfers 1 and 3. (B) N2 versus roam-1MY14 NIL (C) N2 srx-43(lf) versus roam-1MY14 srx-43(lf) (D) N2 daf-22(lf) versus roam-1MY14 daf-22(lf), without or with exogenous icas#9 (10 nM).

E) “Patchy lawn” competition between N2 and roam-1MY14 NIL.

Grey points = individual competition experiments, red line = mean, ***P<0.001, **P<0.01, *P<0.05 compared to expected value of 0.5 by t-test with Bonferroni correction; ns, not significant.

The two tested strains differ in the 182 kb roam-1 region encompassing 81 genes. To ask if the competitive advantage required srx-43, the experiment was repeated using N2 srx-43(lf) and roam-1MY14 srx-43(lf) strains. The competitive N2 advantage disappeared in this setting, identifying the icas#9 receptor SRX-43 as essential for the competitive effect (Figure 5C).

The role of endogenous pheromones was assessed by repeating the competition experiments with N2 and roam-1MY14 strains mutant for the gene daf-22, which is required for the secretion of ascarosides including icas#9. daf-22 mutations eliminated the competitive advantage of the N2 strain over roam-1MY14 (Figure 5D). A competitive advantage of N2 daf-22 was partially recovered in the presence of exogenous 10 nM icas#9 (Figure 5D). These results show that selection on the roam-1 locus depends on pheromones.

The increased roaming of roam-1MY14 at high density might be expected to cause greater exploration of a patchy food environment. In a second competition design, 20 N2 and 20 roam-1MY14 adults were used to seed a patchy environment consisting of 16 small bacterial lawns (Figure 5A). In these conditions the N2 advantage was lost, and instead a moderate but significant selection favored roam-1MY14 over N2 animals (Figure 5E). Together, these results demonstrate that roam-1 can affect fitness bidirectionally depending on srx-43, pheromone production, and food distribution.

DISCUSSION

Conspecific individuals are informative elements of an animal's natural environment, in part because they compete for resources. Our results demonstrate that conspecific pheromones alter long-term foraging strategies, and that natural variation in this behavior stems from altered expression of the icas#9 receptor SRX-43. The complement of ascarosides produced by C. elegans varies with sex, age, and feeding status[27,28], and the specificity of receptors such as SRX-43 provides a mechanism by which this information can be detected by the nervous system to regulate different behaviors and physiological responses. srx-43 is expressed in ASI sensory neurons, which are also targets of internal neuromodulators that regulate roaming and dwelling[8], providing a site of integration of internal and external influences on foraging behavior. Although we do not know the suite of pheromones that are produced by C. elegans in the wild, the presence of secreted icas#9 in dense culture supernatants at concentrations 100-fold above those that suppress roaming suggest that it is a relevant regulator of foraging, and that altered sensitivity to this molecule could affect animals’ overall sensitivity to secreted ascarosides.

The roam-1 QTL that encompasses srx-43 has sequence features of an area under balancing selection. While srx-43 is an essential component of roam-1, it may not be the only gene in this QTL or the only gene under balancing selection, as the haplotype extends for ~30 kb to include several other genes. Moreover, the behavior identified here need not be the most important one in natural settings; it may represent one of several behavioral and physiological responses that facilitate adaptation to different environments. Balancing selection may be fairly common throughout the C. elegans genome: a recent report identified 61 highly divergent regions likely to be under balancing selection that segregate among wild strains of C. elegans, including a second region 200 kb from roam-1[24,29](Extended Data Figure 8). The composition of these regions is biased toward particular gene classes including chemoreceptors, which may act as hotspots of evolution.

Extended Data Figure 8

Recombination between srx-43 and glc-1 in natural isolates

a) The glc-1 gene has previously been shown to be under balancing selection[33], and is near srx-43. The blue line shows SNPs/kb for N2 and MY14 averaged over 5 kb intervals for the region spanning srx-43 and glc-1. The large region of low heterozygosity between srx-43 and glc-1 indicates that balancing selection on glc-1 is unlikely to account for the high heterozygosity near srx-43.

b) Dendrogram for the glc-1 region for strains shown in Figure 4a. roam-1MY14 clades and glc-1 clades are not identical.

At a conceptual level, behavioral genetics in animals including humans is dominated by evidence for gene-environment interactions[30,31]. Our results take this abstraction to a concrete level, showing that natural trait variation acts explicitly at the intersection of innate circuits and environment cues, with genetic changes allowing differential incorporation of environmental information into innate foraging behaviors.

METHODS

Nematode Culture

All strains were grown at 21-22°C on nematode growth medium plates seeded with Escherichia coli OP50 bacteria[32]. For OP50 cultures a single colony was inoculated into 100 ml of LB and grown for 48 hours at 21-22°C. Transgenic lines were generated by standard injection methods, and included the desired transgene, a fluorescent co-injection marker, and empty vector bringing the total DNA concentration up to 100 ng/ul. For each transgene, three independent extrachromosomal lines that propagated the transgene at high rates were tested in parallel to account for variability typical of such strains. All mutagenized strains were backcrossed 5-7 times before characterization.

Strains

MY14-CX12311 RILs

CX14697-CX14712, CX14731-CX14748, CX14750-CX14757, CX14783, CX14784, CX14786-CX14820, CX14822-CX14839. Genotypes inferred from low-coverage genomic sequence and behavioral data are included as Supplementary Table 1.

Near-Isogenic Lines

CX15881 kyIR142 [V:~14.3-~16.8Mb, CX12311>MY14]

CX15878 kyIR139 [V:~14.3-~16.8Mb, MY14>CX12311]

CX15883 kyIR144 [V:~14.3-~16.8Mb, MY14>N2]

CX16075 kyIR147 [V:~15.861-~16.8Mb, MY14>N2]

CX16140 kyIR153 [V:~16.043-~16.8Mb, MY14>N2]

CX16300 kyIR163 [V:~15.861-~16.043Mb, MY14>N2] “roam-1MY14”

CX16294 kyIR157 [V:~15.861-~16.006Mb, MY14>N2]

Transgenic Lines

CX16884 kyIR163 V; kyEx5851 [Psrx-43::srx-43::sl2::GFP @ 2.5 ng/ul, Pmyo3::mcherry @ 5 ng/ul]

CX17202 kyIR163 V; kyEx6012 [Psrx-43::srx-43(nonsense)::sl2::GFP @ 2.5 ng/ul, Pmyo3:: mcherry @ 5 ng/ul]

CX16881 srx-43(gk922634) V; kyEx5848 [srx-43 @ 2.5 ng/ul, Pmyo3::mcherry @ 5 ng/ul]; gk922634 changes R160 to an opal stop codon.

CX17204 kyEx6013 [Psrx-43(N2)::srx-43(N2)-GFP translational fusion @ 50 ng/ul, Pelt-2::GFP @ 5 ng/ul]

CX16943 kyIR163 V; kyEx5894 [Psrx-43(MY14)::srx-43(MY14)::sl2::GFP @ 2.5 ng/ul, Pmyo3:: mcherry @ 5 ng/ul]

CX16425 kyIs602 [Psra-6::GCaMP3.0 @ 75ng/uL; Pcoel::GFP @ 10 ng/uL]; kyEx5594 [Psra-6::srx-43(N2) @ 50 ng/ul, Pmyo3::mcherry @ 5 ng/ul]

CX16931 kyIs602; kyEx5885 [Psra-6::srx-43(MY14) @ 50 ng/ul, Pmyo3::mcherry @ 5 ng/ul]

CX17196 kySi66 [MosSCI Psrx-43(N2)::srx-43(N2)] II; srx-43(gk922634) V, outcrossed 4X

CX17198 kySi68 [MosSCI Psrx-43(MY14)::srx-43(MY14)] II; srx-43(gk922634) V, outcrossed 4X

CX17201 kySi71 [MosSCI Psrx-43(N2)::srx-43(MY14)] II; srx-43(gk922634) V, outcrossed 4X

CX17203 kySi72 [MosSCI Psrx-43(MY14)::srx-43(N2)] II; srx-43(gk922634) V, outcrossed 4X

FK181 ksIs2 [Pdaf-7::GFP + rol-6(su1006)]

CX16958 kyIR163 V; ksIs2

Mutants

CX16849 srx-43(gk922634) V outcrossed 5X to N2. gk922634 is R160opal. This mutation was generously provided by the Million Mutation project[23].

CX16935 kyIR163 srx-43(ky1019) V. ky1019 is a Crispr/Cas9-induced indel mutation that causes a frameshift mutation after the first transmembrane domain (insertion (tcactgagttcgaat), deletion (CCCCG), final sequence TCGCAGCTCTCAAGTtcactgagttcgaatTTCGGAATTCTC). We used the coCRISPR protocol developed in Arribere et al 2014[33]. Young adults were injected with a mix of plasmids containing Cas9, guideRNA targeting rol-6, and guideRNA targeting the location of the desired mutation, as well as a ssDNA template for inducing a dominant rol-6(su1006) mutation. F1 animals with a roller phenotype were isolated and allowed to lay eggs before secondary screening for the target mutation by Sanger sequencing.

JT5464 daf-7(e1372) III; daf-3(e1376) X

CX17307 daf-7(e1372) III; kyIR163 V; daf-3(e1376) X

CX13846 daf-22(ok693) II

CX17082 daf-22(ok693) II; kyIR163 V

Behavioral Analysis

Exploration assays[8] were conducted on 35 mm Petri dishes evenly seeded with 100 μl of OP50 bacteria 24 hours before the start of the assay. Individual two-day old L4 hermaphrodites were picked to the center of the plate. After 16 hours, plates were placed on a grid containing 35 mm squares, and the number of full or partial squares containing tracks were quantified by an investigator blind to the genotype. Pheromones or control solvent were mixed into the agar. A pheromone response for each animal on an ascaroside plate was determined with respect to the behavior of control animals that were tested on ascaroside-free plates on the same day.

Individual pheromone response = (mean # squares entered by controls tested on the same day) – (# squares entered by test animal).

Group pheromone response = mean pheromone response of all individuals tested across days. For statistical analysis, n = total number of animals in the group.

N2-derived strains were tested in 21% oxygen (Figures 1,2dii,iii, 3,5. Extended Data Figures 1,5). All naturally isolated strains and CX12311-derived strains bearing ancestral alleles of npr-1 and glb-5 were tested in 8% oxygen to suppress the oxygen-dependent roaming behavior of ancestral npr-1 alleles[12,34] (Figures 2a,b,di,4, Extended Data Figure 2),

Direct examination of roaming and dwelling was modified from Flavell et al 2013[8]. 14.5 hours before testing, 25 L4 larvae were picked to 150 mm test plates thinly seeded with 1.5 mL of OP50 bacteria with or without synthetic pheromone. Video recording was conducted under red light to minimize behavioral response to imaging conditions. 1.5 hour-long videos were recorded at 3 frames/s using Streampix software (Norpix Inc., Montreal, CA) and a 6.6MP PL-B781F CMOS camera (PixeLINK, Ottawa, CA). Custom Matlab scripts[8] were employed to determine worm trajectories and conduct a two-state hidden Markov model determining the most probable state path for each animal and thereby measure roaming and dwelling state durations.

The low basal exploration rate in daf-7(lf) mutants[7] prevented a direct assessment of the effect of icas#9 on foraging behavior. Instead, we examined daf-7 daf-3 double mutants, as daf-7 canonically acts by antagonizing daf-3, which encodes a co-SMAD transcriptional regulator. daf-3 mutations suppressed the low basal exploration rate of daf-7 mutants. N2 daf-7(lf) daf-3(lf) animals explored control plates moderately more than wild-type, so larger (10 cm) exploration assay plates were used to score these strains.

Statistics

Most experiments were repeated on three separate days. For exploration assays, the standard group size on a single day was 6; this ensured sufficient power to detect moderate effects, while also limiting the influence of daily variation. All plates with a healthy adult animal at the end of the assay were scored and included in the analysis. Randomization was ensured using the following or similar approach: at the start of each exploration assay, 6 animals were placed on a pick at a time, and in the order the animals came off the pick, they were then transferred individually to three control plates and then to three icas#9 plates. Assays were scored by an experimenter blind to the condition or genotype.

Most statistical comparisons were done by ANOVA with Dunnett correction for multiple comparisons or (two-sided) t-test, as noted in the figure legends. The normality of the data was tested with D'Agostino-Pearson omnibus test. Bartlett's test was used to check for differences in variance between groups being statistically compared. N2 groups in 3C and 3K did not pass normality test. As the n was large (>15), ANOVA was still an appropriate test. Moreover, the findings were still significant when a nonparametric test was used (Kruskal-Wallis with Dunn multiple comparison test).

Ascaroside Quantification

150 mL unsynchronized worm cultures were grown for 9 d and fed E. coli (HB101 or OP50), as described[27]. Extracts were generated from the culture medium and analyzed by LC-MS/MS, as described[27], and analyzed on a Thermo Scientific TSQ Quantum Access MAX, with the collision gas pressure set to 1 mTorr. Ascaroside concentrations present in the culture were quantified using the corresponding synthetic standards, except that synthetic ascr#18 was used to quantify ascr#22 and ascr#26, and synthetic icas#3 was used to quantify icas#1 and icas#10.

Recombinant Inbred Lines

The MY14-CX12311 recombinant inbred lines were generated by crossing MY14 males to CX12311 hermaphrodites and CX12311 males to MY14 hermaphrodites, to ensure the mitochondrial DNA from both strains were equally represented in the RILs. 94 F2 were individually picked to plates and inbred through self-fertilization for 10 generations. RIL genotyping was conducted by low-coverage shotgun sequencing[14]. Genomic DNA was fragmented and attached to sequencing adapters with a Nextera DNA Library Prep Kit (Illumina, San Diego, USA). Samples were pooled and sequenced on an Illumina HiSeq 2000. Sequencing reads from each strain were mapped to the WS235 release of the C. elegans genome using bwa to create bamfiles for further analysis[35]. The set of MY14/N2 SNVs identified in the Million Mutation project were used for genotyping purposes[23]. Each genetic variant was genotyped in each strain. Due to the low coverage, the majority of SNVs were not genotyped. To improve the data coverage, we grouped 200 neighboring SNV genotypes together to create a consensus genotype for 540 bins (either N2, MY14 or heterozygous). These genotypes were used for QTL mapping.

QTL Mapping

The pheromone response index was used as the phenotype in combination with the 540 genotype bins from above. R/qtl was used to perform a one-dimensional scan using marker regression on all 540 markers. The significance threshold was determined using 1000 permutation tests. The effect-size of the roam-1 locus was estimated using the fitqtl function with a single QTL. The peak of the roam-1 locus (Chromosome V: 16451686-16579457) was used as an additive and interactive covariate for additional one-dimensional scans, assuming a normal model. The significance threshold for these two tests was also determined using 1000 permutation tests.

NIL Mapping

Before the detailed QTL mapping by sequencing described above, the roam-1 QTL was localized to 2.5 Mb (V:14.3-16.8 Mb) by examining 14 high-confidence phenotypically extreme RILs (Supplemental Table 1). This result, which was confirmed by the full analysis, guided the initial generation of NILs.

kyIR142 was made by backcrossing the RIL CX14816 nine times to MY14, maintaining N2 alleles at V:14.3 and V:16.8 Mb at each generation.

kyIR139 was made by backcrossing the RIL CX14708 nine times to CX12311, maintaining MY14 alleles at V:14.3 and V:16.8 Mb.

kyIR144 was made by crossing kyIR139 with N2 and isolating recombinants with the N2 allele of glb-5 (V:5.56 Mb), the MY14 alleles at V:14.3 and 16.8 Mb, and the N2 allele of npr-1 on X.

kyIR147 and kyIR153 were created by crossing kyIR144 with N2 and identifying progeny with the N2 allele at V:14.3 Mb and the MY14 allele at V:16.8 Mb.

High Density Recombination Mapping

kyIR147 was crossed with males from CX16290, a N2 strain with an integrated fluorescent marker at V:15.83 Mb. F1 progeny were identified by fluorescence, picked to growth plates, and allowed to lay eggs for 12 hours. Following 3 days of growth, ~2600 nonfluorescent F2 were sorted individually into wells of 96-well plates by a worm sorter (COPAS biosort systems; Union Biometrica). These F2 were grown in 200 ul of S-basal buffer + cholesterol supplemented with OP50 bacteria. A fraction of the F3 progeny from each isolate were lysed and genotyped at V:16.043 Mb. Those with an N2 allele at V:16.043 Mb were genotyped at V:15.861 Mb. Twelve recombinants with an N2 allele at V:16.043 Mb and a MY14 allele at V:15.861 Mb were isolated and characterized behaviorally, among which were kyIR163 and kyIR157 (Figure 2Diii). The N2 NIL with kyIR163 (182 kb of MY14 sequence) is referred to as roam-1

Imaging

Calcium imaging experiments were performed and analyzed as described[36]. Briefly, young adult animals were placed into custom-made 3×3 mm microfluidic polydimethylsiloxane devices that permit rapid changes in stimulus solution. Each device contains two arenas, allowing for simultaneous imaging of two genotypes with approximately ten animals each. Animals were transferred to the arenas in S-Basal buffer and paralyzed for 80-100 minutes in 1 mM (−)-tetramisole hydrochloride. Experiments consisted of four 10 s pulses of stimulus separated by 30 seconds of buffer, with 60 additional seconds between stimulus types. Tiff stacks were acquired at 10 frames/second at 5x magnification (Hamamatsu Orca Flash 4 sCMOS), with 10 ms pulsed illumination every 100 ms (Sola, Lumencor; 470/40 nm excitation).

Fluorescence levels were analyzed using a custom ImageJ script that integrates and background-subtracts fluorescence levels of the ASH cell body (4×4 pixel ROI). Using MATLAB, the calcium responses were normalized for each stimulus type by dividing fluorescence levels by the baseline fluorescence, defined as the average fluorescence of the 10 s preceding the first pulse of the stimulus. Each experiment was performed a total of four times over two separate days. Animals of a given strain were pooled together to calculate population mean and standard error (N2 srx-43 allele: 23 animals; MY14 srx-43 allele: 30 animals; array negative control: 19 animals). Experiments were conducted on two days.

For GFP expression studies, live adult animals were mounted on 2% agarose pads containing 5 mM sodium azide. Images were collected with a 100X objective on a Zeiss Axio Imager.Z1 Apotome microscope with a Zeiss AxioCam MRm CCD camera. For daf-7 reporter studies, expression was quantified 16-24 hours after L4 animals were placed on exploration assay plates. Images were processed in Metamorph and ImageJ to generate a maximum intensity Z-projection. Reporter values were assessed as the mean gray value for a 16-pixel-radius circle centered over the cell body minus the mean background intensity. Both ASI neurons were analyzed in each animal; experiments were performed on three days.

Digital PCR

Digital PCR was conducted on a QuantStudio™ 3D digital PCR platform (Thermo Fisher Scientific Inc., New York, USA), and analyzed on the QuantStudio™ 3D AnalysisSuite Cloud. srx-43 mRNA expression studies were conducted on synchronized L4 worms 48 hours after laying. RNA was collected on RNeasy Mini columns (Qiagen, Hilden, DE) and treated with DNAse (Qiagen). SuperScript III First-Strand Synthesis System (Thermo Fisher Scientific) was used to create cDNA libraries. Custom TaqMan Expression Assays (Thermo Fisher Scientific) were used for srx-43 quantification, and the tubulin gene, tbb-1, was used for normalization of digital PCR.

For quantitative analysis of the competition experiments, DNA was extracted with a standard phenol-chloroform protocol. Custom TaqMan SNP Genotyping Assays (Thermo Fisher Scientific) were used to determine the relative ratio of N2 versus roam-1MY14 DNA by digital PCR. The assay was validated with known ratios of N2 to roam-1MY14 DNA (Extended Data Figure 9).

Extended Data Figure 9

Standard curve for dPCR experiments

Best fit line of dPCR results for known ratios of N2 to roam-1MY14 DNA created by mixing different ratios of genomic DNA extracted from independent N2 or roam-1MY14 populations.

Population Genetics

To create the gene and organism phylogenies, we used SNV data downloaded from the Million Mutation project (http://genome.sfu.ca/mmp/) or the CeNDR resource (http://www.elegansvariation.org). For the CeNDR dataset, MY14 was assumed to be clonal or near-clonal with MY23, as was suggested by RAD sequencing. Software was written in Python using the Biopython module to create a neighbor joining tree. For the roam-1 locus, SNVs on Chromosome V between 16,010,000 and 16,030,000 were used. For the glc-1 locus, SNVs on Chromosome V between 16,181,000 and 16,222,000 were used. All SNVs were used to construct the whole genome strain tree. Number of genetic variants and Tajima's D were calculated on 5 kb bins using vcftools[37]. dN/dS was calculated by counting using custom Python scripts analyzing variants between MY23 and the N2 reference. Phylogenies of srx-43 and closely related genes were performed using protein sequences from Thomas and Robertson, 2008[38].

Fitness Assays

Competition experiments consisted of three boom-bust cycles. During the boom phase, population growth led to rapid depletion of food, initiating the bust phase, which lasted for two days. Simple lawn competition experiments were conducted on 100 mm NGM agar plates with a single lawn formed from 800 μl of saturated OP50 culture. Patchy lawn competition experiments were conducted on 150 mm NGM agar plates with a 200 μl ring-shaped OP50 lawn in the center of the plate surrounded by 15 small 40 μl lawns (Figure 5A); at the assay start and at transfers animals were placed in the center of the plate.

Populations were initiated from 20 N2-type and 20 roam-1MY14-type age-synchronized young adult animals. The initial population depleted food within 4 days, and on day 6 animals were washed into M9 media, 20% of the suspension was transferred to a new plate, and the remainder was lysed for quantitative DNA analysis. For the second and third boom bust cycle, food resources were depleted in 2 days and the plates were kept starved for 2 additional days. Following the second bust phase 20% of the animals were transferred to a new plate, and following the third bust phase the entire population was harvested for DNA extraction.

Natural isolates and WT strains	Origin
N2	Bristol, UK
CX12311 kyIR1 [V,CB4856>N2] V; qqIR1 [X, CB4856>N2]	(N2/CB4856)
CX12311 bears ancestral alleles of the npr-1 and glb-5 genes, which affect oxygen sensitivity and are mutated in the N2 laboratory strain[17]; it is therefore used as a comparison strain for wild strains bearing the ancestral alleles.
CB4856 (“HW”)	Hawaii, USA
JU258	Riberio Frio, Madeira
MY1	Lingen, Germany
MY14	Mecklenbeck, Germany
JU775	Lisbon, Portugal
JU1400	Sevilla, Spain
JU1652	Montevideo, Uruguay
AB1	Adelaide, Australia
MY16	Mecklenbeck, Germany
JU1171	Concepcion, Chile
MY6	Roxel, Germany
JU360	Franconville, France
ED3021	Edinburgh, Scotland
MY2	Roxel, Germany