Linked supergenes underlie split sex ratio and social organization in an ant.

Sexually reproducing organisms usually invest equally in male and female offspring. Deviations from this pattern have led researchers to new discoveries in the study of parent-offspring conflict, genomic conflict, and cooperative breeding. Some social insect species exhibit the unusual population-level pattern of split sex ratio, wherein some colonies specialize in the production of future queens and others specialize in the production of males. Theoretical work predicted that worker control of sex ratio and variation in relatedness asymmetry among colonies would cause each colony to specialize in the production of one sex. While some empirical tests supported theoretical predictions, others deviated from them, leaving many questions about how split sex ratio emerges. One factor yet to be investigated is whether colony sex ratio may be influenced by the genotypes of queens or workers. Here, we sequence the genomes of 138 Formica glacialis workers from 34 male-producing and 34 gyne-producing colonies to determine whether split sex ratio is under genetic control. We identify a supergene spanning 5.5 Mbp that is closely associated with sex allocation in this system. Strikingly, this supergene is adjacent to another supergene spanning 5 Mbp that is associated with variation in colony queen number. We identify a similar pattern in a second related species, Formica podzolica. The discovery that split sex ratio is determined, at least in part, by a supergene in two species opens future research on the evolutionary drivers of split sex ratio.

The relative investment in male versus female offspring is a vital fitness component of sexually reproducing organisms. Research on sex allocation theory has yielded breakthroughs in our understanding of topics as diverse as parent–offspring conflict, evolution of cooperative breeding, and genomic conflict (1).

Among these three topics, parent–offspring conflict is predicted to occur in subdivided populations with strong local mate competition, as seen in polyembryonic parasitoids (2), and in systems with relatedness asymmetry between sisters and brothers, as found in haplodiploid species such as the primitively eusocial wasp Polistes chinensis antennalis (3). Considering the evolution of cooperation, parental control of sex ratio is thought to contribute to the maintenance of cooperative breeding; for example, Seychelles warblers living in high quality territories where helpers provide strong benefits produce an excess of females, the helping sex (4). However, similar patterns of biased sex allocation increasing the frequency of the helping sex are not found among all cooperatively breeding birds (5). The first clear empirical example of intragenomic conflict was based upon the discovery of a chromosome that skews sex ratio from female biased to 100% male in the jewel wasp Nasonia vitripennis (6). This paternal sex ratio chromosome is transmitted through sperm to fertilized eggs, where it causes the loss of other paternally inherited chromosomes to produce exclusively male offspring (7, 8). Subsequent discoveries of sex ratio distorter systems unfolded in different directions, including female-biased sex ratios mediated by endosymbionts (9, 10). These studies opened the door for additional research on intragenomic conflict in multiple contexts, including between sexes (11, 12) and between social insect castes (13).

Where there is intragenomic conflict, one resolution is evolution of suppressed recombination to reduce the frequency of deleterious multilocus genotypes. This is illustrated in the standard model of sex chromosome evolution (14, 15), in which selection favors the loss of recombination between a sexually antagonistic locus and a sex-determining locus on the same chromosome, eventually leading to a Y or W chromosome that is exclusively present in one sex. Under the “reduction principle” (16), this is also expected to occur around sex-ratio distorters. In line with this prediction, sex-ratio distorter loci often occur in regions of low recombination (17–20), but we lack evidence for the direction of causality. The reduction principle is also expected to contribute to the formation of autosomal supergenes controlling other complex traits that involve epistatic interactions between two or more loci. Such supergenes have been found to control phenotypes including polymorphic wing coloration in butterflies (21), mating strategies in birds and fungi (22–25), self-incompatibility in plants (26), and colony social organization in ants (27, 28). Autosomal supergenes, like sex chromosomes, are likely to represent the resolution of past intragenomic conflict between two or more loci.

Supergenes underlie at least two independently evolved cases of social polymorphism in ants. In the fire ant Solenopsis invicta, colony queen number is controlled by a supergene spanning most of a single chromosome (27). Formica selysi has a similar chromosome-spanning supergene underlying colony queen number, but there is no detectable overlap in gene content between the two (28). More recently, both ant social supergenes were shown to underlie colony queen number in other congeneric species (29, 30). In both systems, the haplotype associated with multiqueen (= polygyne) social structure is a selfish transmission distorter (31–33). These discoveries raise new questions about links between social structure and sex ratio that have been proposed in classic literature about sex allocation in Hymenoptera.

Trivers and Hare (34) proposed that queen–worker conflict, which is shaped by relatedness asymmetry within each nest, drives biased sex ratios. Sex ratios represent the proportion of reproductive females and males and do not include workers. Since workers are more related to their full sisters (average relatedness = 0.75) than to their brothers (average relatedness = 0.25), workers in single-queen, monandrous colonies should favor the production of queens over males. Trivers and Hare (34) predicted that worker interests would prevail in these cases since workers provide all care for the brood and that this would result in female-biased offspring production. Queens are equally related to male and female offspring, so they should generally favor a 1:1 sex ratio. In colonies with multiple queens or a single, multiply inseminated queen, the lower relatedness reduces this conflict between queens and workers, resulting in weaker selection for biased sex allocation (1, 34). Although these predictions revolutionized the way that researchers think about how relatedness shapes inclusive fitness in social insect colonies, they are not ubiquitously upheld in empirical studies (1).

Strikingly, some social insect species exhibit a nearly complete segregation of male and queen production at the colony level, in a phenomenon known as “split sex ratio.” Such extreme cases have been observed in at least 20 different genera of ants, wasps, and bees (35–37). Boomsma and Grafen (36) argued that this pattern is consistent with worker control of sex ratios in ant populations with variation in relatedness asymmetry: workers that are more related than the population average to their nestmates should favor specializing in the production of new queens (hereafter, “gynes”), while those that are less related than average should specialize in male production (36, 38). The variation in relatedness asymmetry would normally emerge from the number of mates per queen, or from the number of queens per colony, or both.

The models of Boomsma and Grafen inspired a burst of empirical research on split sex ratios. Subsequent studies (reviewed in refs. 37, 39, and 40) tested four scenarios that increase the variation in relatedness asymmetry among colonies, including variation in queen number (41–46), variation in queen insemination (41, 42, 45, 47), variation in breeder turnover (46, 48–50), and presence or absence of workers (51), as well as two scenarios not involving relatedness asymmetry, namely resource availability (43, 52–54) and maternally inherited parasites (55, 56). Ants in the genus Formica emerged as a prominent model system, as a result of their widespread and well documented variation in sex ratio and social structure (35). Many species exhibit split sex ratios or highly biased sex ratios (34, 41–44, 47) but not all of these examples follow predicted patterns based on relatedness asymmetry. Finnish populations of Formica truncorum and Formica exsecta follow theoretical predictions: in colonies with a single queen (= monogyne), monandrous queens tend to produce gynes, while polyandrous queens tend to produce males (41, 43). A similar pattern was found in monogyne and polygyne colonies in F. truncorum, with polygyne colonies producing males (44). A socially polymorphic population (i.e., comprised of monogyne and polygyne colonies) of F. selysi and a polygynous population of F. exsecta that exhibit variation in relatedness asymmetry deviated from these predicted patterns (45, 46). Additional studies have identified potential roles of habitat and diet in shaping sex allocation in Formica podzolica (43), F. exsecta (53), and Formica aquilonia (54). Other studies have suggested that investment in sexual offspring is mediated by colony needs for queen replacement (48, 57), with gynes being produced by colonies with relatively few queens. Finally, although Wolbachia is present in some Formica species exhibiting split sex ratio, it does not appear to influence sex ratio in any system studied so far (55, 56).

Taken together, it appears that there are yet missing pieces to the puzzle of how and why ants achieve split sex ratios. A meta-analysis attributed about 25% of the observed variance in sex allocation to relatedness asymmetry and variation in queen number (37). Theoretical examinations following from this finding support a possible role for virgin queens (which would produce only male offspring) or queen replacement (58), but another possible factor is that sex allocation by queens is itself under genetic control.

Here, we examine the evidence for genetic control, which could be responsible for some of the unexplained variance in patterns of split sex ratio. We 1) conduct a genome-wide association study (GWAS) for variants associated with sex ratio in Formica glacialis, 2) infer transmission patterns of sex-ratio–associated variants from colony-level genotype frequencies, 3) evaluate whether sex ratio and social organization map to the same region of the genome, and 4) examine the sister species F. podzolica (59) to test for a shared genetic basis of sex ratio.

Results

GWAS of Sex Ratio.

Through a GWAS of 138 F. glacialis whole-genome sequences, we identified numerous variants associated with colony sex allocation in a region of chromosome 3 spanning 5.5 Mbp (Fig. 1 and ). A principal component analysis (PCA) of variants on chromosome 3 revealed three distinct genotype clusters, one of which was observed in just six individuals (Fig. 1). Of the workers with low PC2 scores (yellow and green clusters, Fig. 1), 60.2% were collected from male-producing colonies, while 93.3% of workers with high PC2 scores (purple cluster, Fig. 1) were from gyne-producing colonies. An investigation of genetic differentiation (FST) between genotype clusters on chromosome 3 revealed two adjacent regions of differentiation: between the two clusters with low PC1 scores, we observed differentiation spanning the region from 2 to 7.5 Mbp (Fig. 2), similar to the region revealed in the initial GWAS. Between the two clusters with low PC2 values (both of which harbored an excess of workers from male-producing colonies), we identified a differentiated region from about 7.5 to 12.5 Mbp, as well as a small peak at 2 Mbp (Fig. 2).

Fig. 1.

Alternative haplotypes on chromosome 3 are associated with colony sex ratio in F. glacialis. (A) A GWAS using a linear mixed model implemented in GEMMA reveals a large region on chromosome 3 significantly associated with colony sex ratio in F. glacialis. In this Manhattan plot, each point represents a SNP, with the log-transformed P value from the GWAS plotted against physical position on the genome. (B) A PCA of variants on chromosome 3, with each SNP treated as a variable, identified three clusters corresponding to three supergene genotypes, SmA/SmA (green), SmA/SmD (purple), and Sp/SmA (yellow). In this genotype nomenclature, “S” represents supergene, “m” represents monogyne social structure, “p” represents polygyne social structure, “A” represents male-biased sex ratio, and “D” represents female-biased sex ratio. The shapes represent the observed sex ratio phenotype of the colony each worker belongs to, male-producing (squares), gyne-producing (circles), and undetermined sex ratio (diamond).

Fig. 2.

Genetic differentiation between individuals with alternative genotypes on chromosome 3 reveals two adjacent supergenes in F. glacialis. (A) SmA/SmA workers (green points in Fig. 1) and SmA/SmD workers (purple points in Fig. 1) exhibit elevated FST between 2 Mbp and 7.5 Mbp. (B) SmA/SmA workers and Sp/SmA workers (yellow in Fig. 1) exhibit elevated FST between 7.5 Mbp and 12.5 Mbp, with a small peak at 2 Mbp. Points represent 10-kbp windows.

Association between Social Structure and Supergene.

Previous studies found that colony queen number in F. selysi (28, 60) and other European Formica species (29) is controlled by a social supergene on chromosome 3. To determine whether a supergene on chromosome 3 similarly underlies colony queen number in F. glacialis, we investigated variation in 19 additional colonies from other populations using double-digest restriction-site–associated DNA sequencing (ddRADseq). We calculated opposing homozygosity among nestmates (i.e., the presence of two alternative homozygous genotypes in nestmates) to determine colony social structure (Fig. 3). For each colony, we counted the number of loci that displayed opposing homozygosity within a sample of eight nestmate workers. For any single biallelic SNP, there is no combination of parental genotypes that would produce both of the alternative homozygous genotypes within a set of multiple full-sister offspring (e.g., workers from a monogyne colony). In contrast, it's entirely possible to obtain both homozygous genotypes within a set of workers descended from more than two parents, as expected in a polygyne colony. In practice, monogyne colonies can have low but nonzero levels of apparent opposing homozygosity due to rare genotyping errors or low to moderate levels in colonies headed by a single multiply mated queen. In our dataset, 11 colonies (58%) exhibited very low levels of opposing homozygosity and were inferred to be monogyne, while four colonies (21%) had high levels of opposing homozygosity and were inferred to be polygyne (Fig. 3). Four colonies with intermediate levels of opposing homozygosity were considered to be “ambiguous” for our analyses. The variation in opposing homozygosity mapped to the supergene region (Fig. 3). In particular, single-nucleotide polymorphisms (SNPs) that were significantly associated with variation in social structure in the ddRADseq data were localized in the 7.5 to 12.5 Mbp region (Fig. 3), corresponding to the region identified in Fig. 2.

Fig. 3.

The Sp supergene haplotype is associated with polygyne social structure in F. glacialis from Alaska, British Columbia, and Alberta, based on ddRADseq genotyping of 7 to 8 workers from each of 19 colonies. (A) Opposing homozygosity varies among colonies. Putative monogyne colonies are colored blue, putative polygyne colonies are colored red, and undetermined are colored gray. (B) GWAS reveals multiple SNPs associated with colony-level opposing homozygosity on chromosome 3. (C) Significantly associated SNPs occur within the 7.5 to 12.5 Mbp region, also identified in the Yukon population and shown in Fig. 2. (D) A PCA of variants in both the ddRADseq (filled circles) and whole-genome datasets (open circles) shows that the Sp haplotype identified in the Yukon population (open circles in Upper Right) clusters with the haplotype associated with polygyne social structure in other populations (red). The majority of individuals from the Yukon population (open circles in Left) cluster with workers from monogyne colonies in other populations (blue).

A PCA of markers on chromosome 3 that were shared in the whole-genome and ddRADseq datasets revealed that the colonies that were assessed to be polygyne based on a high frequency of opposing homozygosity (red, Fig. 3) consistently exhibited one genotype. This genotype was shared with the six individuals from the whole-genome sequencing library that formed the yellow cluster in Fig. 1. These individuals were heterozygous for two alternative supergene haplotypes, one of which appears to occur exclusively in polygyne colonies. Based on our ddRADseq results, we inferred that these three colonies (4% of the 71 colonies collected for our split sex ratio analysis) from our focal population are likely to be polygyne. Following the notation developed to describe the monogyne- and polygyne-associated haplotypes in F. selysi (Sm and Sp, respectively, where S represents a supergene and m and p refer to each colony social form; 28), we defined this haplotype as the Sp haplotype of F. glacialis. We note that the Sp found in other Formica species spanned about 10.5 Mbp of chromosome 3, from 2 Mbp to about 12.5 Mbp (29), while the Sp identified here in F. glacialis was shorter. The remaining two genotype clusters identified in the whole-genome dataset (green and purple clusters, Fig. 1) both grouped with workers from colonies assessed to be monogyne in the ddRADseq dataset based on very low levels of opposing homozygosity (Fig. 3). Based on the regions of differentiation among genotype clusters (Fig. 2), we hypothesized that individuals from the purple cluster carried two alternative supergene haplotypes in the 2 to 7.5 Mbp region of chromosome 3 (subsequently confirmed with targeted genotyping; Fig. 4). One of these haplotypes was found almost exclusively in gyne-producing colonies. The other haplotype was usually homozygous in male-producing colonies. Since one genotype was associated with the production of daughters in monogyne colonies, we named these alleles after the mythological twins Danaus and Aegyptus, who respectively had 50 daughters and 50 sons. Individuals from the gyne-producing cluster (purple, Fig. 1) had the genotype SmA/SmD, while those from the predominantly male-producing cluster had the genotype SmA/SmA (green, Fig. 1).

Fig. 4.

Genotype distributions differ between colonies with alternative sex ratios and social structures. (A) The SmA/SmD genotype in F. glacialis occurs in approximately half of workers from gyne-producing colonies and is rare in workers from male-producing colonies. (B) We propose a model of Mendelian inheritance for maintenance of this supergene system in a largely monogyne population. Heterozygous queens (SmA/SmD) mated with an SmA male produce exclusively female offspring with SmA/SmD and SmA/SmA genotypes. Gynes with the heterozygous genotype become gyne producers, while homozygotes become male producers. (C) Gyne-producing colonies usually harbor a mix of SmA/SmA and SmA/SmD workers (31/34 colonies). Male-producing colonies usually contain exclusively SmA/SmA workers and produce SmA males (27/34 colonies). Three additional male-producing colonies contain exclusively Sp/SmA workers and either Sp or SmA males. We did not detect SmA/SmD workers in three gyne-producing colonies, while we found at least one SmA/SmD worker in four male-producing colonies, suggesting that the genetic basis of split sex ratio may be imperfect in this system. (D) Among monogyne colonies from Alaska, British Columbia, and Alberta, two contain exclusively SmA/SmA workers, while nine contain SmA/SmA and SmA/SmD workers. All workers from polygyne colonies carry one Sp haplotype, consistent with the association between the Sp haplotype and polygyne social structure observed in other Formica species (29). The majority of polygyne workers are Sp/SmA, and we detected one Sp/SmD worker.

Distribution of Genotypes in Colonies.

We developed two polymerase chain reaction restriction fragment length polymorphism (PCR-RFLP) assays to distinguish these three genotypes in a larger number of individuals from each of the colonies in the focal population in Yukon Territory. Workers from gyne-producing colonies were a mix of SmA/SmD heterozygotes and SmA/SmA homozygotes, while workers from male-producing colonies were most often SmA/SmA homozygotes or Sp/SmA heterozygotes (Fig. 4). This suggests that gyne-producing monogyne colonies are usually headed by SmA/SmD queens that produce a mix of SmA/SmD and SmA/SmA gynes and workers, while male-producing monogyne colonies are usually headed by SmA/SmA queens (Fig. 4). (Recall that the SmD haplotype is named for Danaus, who had only daughters, while the SmA haplotype is named for Aegyptus, who had only sons.) As a result, the overwhelming majority of males in the population should have the SmA genotype. Looking at each colony, we showed that 31 out of 34 gyne-producing colonies harbored at least 1 SmA/SmD worker out of 8 genotyped, while 27 out of 34 male-producing colonies harbored only SmA/SmA workers and SmA males (Fig. 4). Among the remaining male-producing colonies, three harbored only Sp/SmA workers and either Sp or SmA males (and were likely polygyne). We observed individuals bearing the SmD haplotype in four male-producing colonies. We inferred the genotypes of individuals from colonies with known social structure in the ddRADseq dataset using a set of diagnostic SNPs. Across these additional populations, we showed that two monogyne colonies harbored exclusively SmA/SmA workers, while nine harbored a mix of SmA/SmD and SmA/SmA workers (Fig. 4). The four polygyne colonies all contained Sp/SmA workers; one colony contained a single Sp/SmD worker as well.

Homologous Supergene in Sister Species.

We obtained a smaller sample of colonies of F. podzolica, the sister species of F. glacialis (59), that exhibited split sex ratios at the focal site in the Yukon Territory. While the GWAS analysis was inconclusive (), we observed similar qualitative patterns in the genomic differentiation between genotype clusters identified in a PCA (Fig. 5). Individuals from these two PCA clusters (Fig. 5) exhibited elevated genetic differentiation from 2 to 7.5 Mbp along chromosome 3 (Fig. 5). Gyne-producing colonies harbored a mix of putative SmA/SmD heterozygotes and SmA/SmA homozygotes. The majority of male-producing colonies contained exclusively SmA/SmA workers (Fig. 5). A large number of SNPs distinguishing SmA and SmD haplotypes were conserved between F. podzolica and F. glacialis (Fig. 5).

Fig. 5.

Alternative genotypes matching those detected in F. glacialis are also associated with colony sex ratio in Formica podzolica. (A) A PCA of genetic markers on chromosome 3 reveals two genotype clusters, SmA/SmA (green) and SmA/SmD (purple). (B) In F. podzolica, the region of elevated FST between SmA/SmA and SmA/SmD workers spanned 2 to 7.5 Mbp along chromosome 3. Points represent 10-kbp windows. (C) Among workers from gyne-producing colonies, about half carry the SmA/SmA genotype and half carry the SmA/SmD genotype. In contrast, workers from male-producing colonies are almost exclusively SmA/SmA. (D) The F. glacialis and F. podzolica SmD haplotypes associated with female-biased sex ratio appear to be homologous, based on the high concentration of SmD-specific SNPs shared between the two species in the 2 to 7.5 Mbp region of chromosome 3. Each vertical bar in this panel represents a 10-kbp window on chromosome 3, with the height of the bar denoting the number of SmD-specific SNPs shared between the two species in that window.

Discussion

We demonstrate that a chromosome underlying queen number across the Formica genus is also associated with the split sex ratios observed in a sister species pair. Sex ratio variation based on queen genotype could account for some of the empirical exceptions (reviewed by ref. 37) to the patterns predicted by Boomsma and Grafen (36, 38). In F. glacialis, we show that the SmD supergene haplotype behaves like a “W” sex chromosome in that it is present almost exclusively in females and in a heterozygous state. A key difference is that it influences the sex ratio of reproductive offspring rather than the sex of the individual bearer, as does the W’ chromosome of the Hessian fly (61). Single-queen gyne-producing colonies generally harbor a mix of SmA/SmD and SmA/SmA workers, suggesting that the queens are SmA/SmD heterozygotes crossed with SmA males. Through Mendelian inheritance, half of their reproductive daughters (the heterozygotes) will in turn be gyne-producing queens, while the other half will be male producers. Males are produced either by homozygous SmA/SmA single queens or by polygyne (Sp/SmA) queens. We noted a few exceptions to this pattern in both gyne- and male-producing colonies. These exceptions could indicate that genetic control is imperfect, that some colonies that truly produced both males and gynes were misclassified due to our sampling effort taking place on a single day, or they could result from other factors affecting the field colonies.

Linked Supergenes Underlie Sex Ratio and Queen Number.

A striking finding of this study is that the overall length of the social supergene discovered in other Formica species appears to be split into two adjacent, linked supergene regions in F. glacialis. One half of the supergene, from 2 to 7.5 Mbp on chromosome 3, is associated with split sex ratio. The other half, from 7.5 to 12.5 Mbp, which includes the gene knockout identified as a candidate conserved gene influencing social structure in other Formica species (29), is associated with social structure (Fig. 2).

Possible Evolutionary Origins of Linked Sex Ratio and Social Supergenes.

Theory predicts split sex ratio to evolve in social hymenopteran populations with variation in relatedness asymmetry. Starting from this point, we propose two possible scenarios that could explain the evolution of these linked regions. In one scenario, we speculate that split sex ratio may have evolved in socially polymorphic Formica populations, wherein monogyne and monandrous queens would specialize in gyne production, while polygyne or polyandrous queens would produce predominantly males. Such patterns were documented in other Formica species, including F. truncorum (47, 62) and Finnish populations of F. exsecta (63), although we note that this pattern is not present in all previously studied Formica species (45, 48). Specialization in offspring sex ratio based on social structure would select for reduced recombination between loci influencing sex ratio and social structure. In populations with little relatedness asymmetry, as observed in our predominantly monogyne F. glacialis population in the Yukon, rare recombinant supergene haplotypes that decouple social determination from sex ratio determination could spread in the population. In this case in particular, we suggest that recombination or gene conversion may have resulted in the transfer of supergene regions that evolved on the Sp haplotype to the Sm background, leading to the formation of the SmA haplotype associated with monogyne social structure and the production of males. Such sex ratio supergene systems may persist in species with a mix of socially polymorphic and socially monomorphic populations. This genetic control could explain deviations from the theoretical predictions of Boomsma and Grafen (36). Deviations from theoretical predictions have been found in both socially polymorphic and socially monomorphic populations: in a socially polymorphic F. selysi population, one social form exhibits strongly split sex ratios and the other is intermediate (45), while in a polygyne F. exsecta population, most colonies produce an excess of male offspring (46).

In an alternative scenario, we speculate that a gene or supergene influencing sex ratio could predate the appearance of persistent social polymorphism; when alternative social structures emerged, selection for male-biased production in colonies with lower average relatedness and for gyne-biased production in colonies with higher average relatedness could have led to the appearance of linked genetic variants favoring one or more queens. The dual roles of linked supergenes in shaping social organization and sex ratio in Formica species could help to explain why this supergene has persisted for millions of years (29). Future studies could examine these speculative scenarios by seeking evidence of sex ratio supergenes in other, distantly related Formica species.

Sex-Specific Genetic Variation in Ants.

Our study in F. glacialis is not the first to identify sex-specific genetic differences between ant gynes and males (20, 64). However, the mechanisms that produce these sex-specific genetic differences appear to differ across systems, and none are fully understood. Kulmuni and Pamilo (64) showed that hybridization between F. aquilonia and Formica polyctena results in admixed females but that surviving males tend to have a genotype comprised of alleles from only one parental species (64, 65). They proposed that recessive incompatibilities between the genomes of the two species are exposed to selection in haploid males (but see ref. 66). In the tawny crazy ant Nylanderia fulva, males invariably carry the same allele at 2 out of 12 microsatellite loci, while females are almost always heterozygous at these loci (20). Diploid N. fulva eggs that were homozygous for the male-associated alleles at these loci failed to develop. In F. glacialis, we similarly observe a haplotype, SmD, that is found almost exclusively in females. However, in the aforementioned cases, certain genotypes have sex-specific effects on viability, while in the case of F. glacialis and F. podzolica, genotypes affect which sex is produced (proposed inheritance outlined in Fig. 4). In this way, sex-specific lethality of certain genotypes is not necessary to explain the observed pattern in F. glacialis, although we cannot rule out such mechanisms. Further research is needed to identify the mechanisms maintaining these genetic differences between the sexes in all systems.

Reconciling Environmental and Genetic Influences on Sex Ratio.

Some previous empirical discoveries still need further examination in the light of the newly discovered split sex ratio supergene. Several experimental studies provided evidence that environmental quality and diet can influence colony sex ratio, including in a population of F. podzolica from central Alberta (43). Given the strong genetic effect that we uncovered, the relative roles of the environment and genetics should be reconciled in future studies. Although the linkage between social and sex ratio supergenes hints at a role of parent–offspring conflict in shaping split sex ratio in Formica ancestors (34, 38), many questions remain about how worker control could function in a system with genetic determination of sex ratio. Understanding how the sex ratio supergene functions will help to illuminate how the contemporary conflict plays out. For example, does the SmD haplotype cause cessation of haploid egg production? What factors prevent female offspring of SmA/SmA queens from developing into gynes instead of workers? Can workers override genetic control by manipulating larval nutrition?

Comparison of Formica Social Supergenes.

Despite limited sampling, we also document an intriguing deviation in the mode of action of the F. glacialis social supergene compared to that of F. selysi. Across the individuals in the ddRADseq dataset collected from polygyne colonies, we did not detect any Sp/Sp homozygous individuals (Fig. 4). Our sample of polygynous colonies harbored almost exclusively Sp/SmA workers (n = 30), with only a single Sp/SmD worker. In contrast, previous studies of F. selysi found that polygynous colonies contained exclusively Sp/Sm and Sp/Sp workers and Sp males (28, 33, 60). These studies detected no systematic variation at the supergene in monogynous F. selysi workers (28, 60), but we note that analyses were carried out with relatively sparse ddRADseq markers, so it is possible that a sex ratio supergene could be present in F. selysi.

Conclusion.

Here, we describe a supergene that appears to have a significant influence on offspring sex ratio in ants. This sex ratio supergene is closely linked with a previously described supergene that underlies colony queen number in Formica ants. The discovery that split sex ratio can have a genetic basis may help to resolve the conflicting empirical results about whether and how split sex ratio emerges to resolve parent–offspring conflict in social hymenopterans. We suggest that genetic control of sex ratio should be investigated in other social insects, particularly in those that do not conform to theoretical predictions.

Materials and Methods

Sampling and Field Observations.

We sampled a mixed population of F. glacialis and F. podzolica 50 km west of Whitehorse, Yukon Territory, Canada on July 17, 2016. We removed the top 5 to 10 cm of soil from each nest mound and assessed the presence and sex of winged sexuals. Approximately 50% of colonies examined had sexuals present. When we observed strongly biased sex ratios (i.e., of the first 10 sexuals examined, at least 9 were of the same sex), we sampled at least 8 workers and up to 5 males. We did not sample gynes. In total, we sampled 71 F. glacialis colonies, of which 34 were male producing and 34 gyne producing. The remaining three sampled colonies contained no F. glacialis sexuals; two of the three also contained workers of the socially parasitic species Formica aserva. We estimate that 80 to 90% of colonies with winged sexuals exhibited a biased sex ratio. Because we observed a sample of sexuals at a single point in time, the true lifetime sex ratio of some colonies may differ from the sex ratio we observed.

Whole-Genome Sequencing.

We sequenced 138 genomes of workers from 71 F. glacialis colonies. We extracted genomic DNA using the Qiagen DNeasy insect tissue protocol and prepared whole-genome DNA libraries using a low-volume Illumina Nextera protocol (67) as modified by ref. 68. The libraries were sequenced on an Illumina HiSeq X-Ten by Novogene, Inc., using 150 base-pair (bp) paired-end reads.

Variant Calling.

We merged overlapping paired-end reads with PEAR version 0.9.10 (69), aligned the reads to the F. selysi reference genome (29) using BWA-MEM version 0.7.17 (70), and removed PCR duplicates with Samtools version 1.8 (71). We called variants using Samtools mpileup version 1.8 (72) and filtered the genotypes for missing data (20% per locus, −max-missing 0.8), minor allele count (−mac 2), and minimum depth (−minDP 1) with VCFtools version 0.1.13 (73).

Population Genetic Analyses.

We identified regions significantly associated with colony sex ratio in 71 F. glacialis colonies (n = 138 individuals) by performing a GWAS using a univariate linear mixed model implemented in GEMMA version 0.94 (74), with significance assessed using a Wald test. We adjusted P values to account for multiple comparisons using the false discovery rate (FDR), with an FDR-adjusted significance threshold of 0.05. This model uses a genetic similarity matrix to control for population structure as a random effect. Colony sex ratio was coded as a binary phenotype, with workers from male-producing colonies coded 0 and from gyne-producing colonies coded 1. Workers from three colonies without a sex ratio phenotype were assigned an “NA” phenotype. We visualized these results with Manhattan and quantile–quantile plots using the qqman package in R (75). Based on the GWAS results showing a large region of chromosome 3 significantly associated with sex ratio, we performed a PCA using Plink version 1.90b3.38 (76) on variants on this chromosome, which resulted in three distinct genotypic clusters that we named SmA/SmA, SmA/SmD, and SmA/Sp (see Results). In order to identify the regions of the chromosome that are responsible for the differences between these clusters, we calculated Weir and Cockerham's FST (77) in 10-kbp windows between the three clusters using VCFtools version 0.1.13 (67).

Comparisons with Sister Species F. podzolica.

We examined the underlying genetics of split sex ratio in the sister species F. podzolica, as well. We sampled 12 colonies (5 male-producing, 7 gyne-producing) from the same Yukon locality and sequenced the genomes of 22 workers and called variants using the same methods as for F. glacialis. We identified genetic clusters using a PCA based on variants on chromosome 3 and calculated FST between genetic clusters, again using the methods described in Population Genetic Analyses for F. glacialis. In order to assess homology between the alternative supergene haplotypes of the two species, we identified SNPs with alleles specific to the SmD haplotypes of both species by comparing allele frequencies across chromosome 3 in four groups: F. glacialis SmA/SmA, F. glacialis SmA/SmD, F. podzolica SmA/SmA, and F. podzolica SmA/SmD. Loci with putative SmD-specific alleles shared between both species were defined as those which have allele frequency between 0.4 and 0.6 in both of the SmA/SmD groups and allele frequency >0.95 or <0.05 in both of the SmA/SmA groups. We plotted the frequency of SNPs meeting these criteria in 10-kbp windows along chromosome 3; regions containing a high frequency of transspecies SmD-specific SNPs provide evidence for homology of the SmD haplotypes in the two species.

Assessment of Colony Social Organization.

We sampled 8 workers from 19 additional F. glacialis colonies in Alaska, British Columbia, and Alberta, where no winged sexuals were visible at the time of collection. We genotyped 145 of these workers using the ddRADseq protocol of Brelsford et al. (78), with restriction enzymes SbfI and MseI. These libraries were sequenced on the Illumina HiSeq 4000 platform by the QB3 Genomics core facility of University of California Berkeley, with 100-bp paired-end reads. We aligned reads and called variants using the procedures described in Variant Calling for whole-genome data but omitting the removal of PCR duplicates. Raw variants were filtered using VCFtools version 0.1.13 (73), retaining genotypes with sequence depth of at least 7 and variants with genotype calls in at least 80% of samples.

To assess social organization of these 19 colonies, we calculated the number of loci exhibiting opposing homozygosity within each colony: that is, at least one worker homozygous for the reference allele and one worker homozygous for the alternate allele. In haplodiploid organisms, a male transmits the same allele to all of his offspring, so in a group of full siblings, opposing homozygosity is expected to be absent except in the cases of genotyping errors or de novo mutations. Colonies with multiple queens, or with a multiply mated queen, are expected to have a higher number of loci with opposing homozygosity.

We conducted a GWAS for variants associated with colony-level opposing homozygosity using a linear mixed model implemented in GEMMA version 0.94 (74), which uses a relatedness matrix to control for nonindependence of samples. Finally, we carried out a PCA of variants on chromosome 3 on a merged dataset of whole-genome and ddRAD F. glacialis genotypes. We generated a list of variants on chromosome 3 present in both the whole-genome and ddRADseq filtered variant call format (VCF) files, extracted those variants from both datasets, and generated a merged VCF using VCFtools version 0.1.13 (73). We used Plink version 1.90b3.38 (76) to carry out a PCA on the resulting merged genotypes. The PCA revealed clusters of individuals with the same genotype across the merged whole-genome and ddRAD datasets.

Colony Genotype Distributions.

We designed a targeted PCR-RFLP assay for a transspecies SNP tagging alternative chromosome 3 haplotypes in both F. glacialis and F. podzolica. We designed primers (CTGGAACAACGGATCCTCA and TTCGCGATTCGAATTTCTC) to amplify a 338-bp fragment, which, when digested with the restriction enzyme MluCI, produces fragments of 325 and 13 bp for the haplotype associated with gyne production and 223, 102, and 13 bp for the haplotype associated with male production. We used this assay to genotype six additional workers and any available males for all colonies of both species.

We designed a second PCR-RFLP assay for a transspecies SNP broadly conserved across Formica that differs between Sm and Sp alleles of the gene knockout. Primers GGTGGYTCTTTCAACGACG and GCCATGTTCACCTCCACCA amplify a 230-bp fragment, which when digested with the restriction enzyme HinfI produces fragments of 132 and 98 bp for the Sm allele and 230 for the Sp allele. We used this assay to genotype six additional workers and any available males from three F. glacialis colonies where initial whole-genome sequencing of two workers identified the presence of the Sp allele at knockout. For both PCR-RFLP assays, we visualized the distinct banding patterns with 2% agarose gel electrophoresis.

We constructed bar plots of supergene genotype frequency by colony based on whole-genome sequencing and PCR-RFLP genotyping for the Yukon population and based on ddRAD for Alaska, Alberta, and British Columbia populations.