Determining social and population structures requires multiple approaches: A case study of the desert ant Cataglyphis israelensis

Abstract The remarkable diversity of ant social organization is reflected in both their life history and population kin structure. Different species demonstrate a high variation with respect to both social structure and mating strategies: from the ancestral colony type that is composed of a single queen (monogyny), singly inseminated (monoandry), to the more derived states of colonies headed by a multiply inseminated queen (polyandry), to colonies composed of multiple queens (polygyny) that are either singly or multiply inseminated. Moreover, the population structure of an ant species can range from multicoloniality to polydomy to supercoloniality, and Cataglyphis is considered to be a model genus in regard to such diversity. The present study sought to determine the social and population structure of the recently described C. israelensis species in Israel. For this purpose we employed a multidisciplinary approach, rather than the commonly used single approach that is mostly based on genetics. Our study encompassed behavior (nest insularity/openness), chemistry (composition of nestmate recognition signals and cuticular hydrocarbons), and genetics (microsatellite polymorphism). Each approach has been shown to possess both advantages and disadvantages, depending on the studied species. Our findings reveal that C. israelensis colonies are headed by a single, multiply inseminated queen and that the population structure is polydomous, with each colony comprising one main nest and several additional satellite nests. Moreover, our findings demonstrate that none of the above‐noted approaches, when employed individually, is suitable or sufficient in itself for delineating population structure, thus emphasizing the importance of using multiple approaches when assessing such complex systems.

, they all share the basic colony structure that is presumably driven by kin selection, reinforced by the haplo-diploid sex determination system, in which within-colony relatedness is higher than that between colonies, and individuals benefit from inclusive fitness rather than individual fitness (Bourke & Franks, 1995;Hamilton, 1964Hamilton, , 1972Hölldobler & Wilson, 1990;Queller, 1992).
However, when there is limited dispersal or variation in the basic colony form, such as multiply mated queens (polyandry) or multiple queens (polygyny), these assumptions do not always hold true.
In such cases, within-colony relatedness is very low (Keller, 1993), and workers from two distinct but neighboring colonies can be more related to each other than to their respective nestmates (Heinze, 2008;Procter et al., 2016;Wright, 1943). It is thus insufficient to assess colony boundaries using genetic methods alone. Another colony-specific "marker" is that of the cuticular hydrocarbon profile (CHC), which in many ant species is responsible for nestmate recognition. CHC composition is attributed both to genetic and nutritional differences (Liang & Silverman, 2000;Vonshak, Dayan, Foucaud, Estoup, & Hefetz, 2009), and nestmates are assumed to possess similar CHC profiles due to the so-called "colony gestalt" (Hefetz, 2007;Soroker, Vienne, Hefetz, & Nowbahari, 1994). However, there may be cases in which these profiles could also overlap with individual profiles of workers from other nests, especially if dispersal is limited. Thus, employing additional methods is important when seeking to understand the colony boundaries. One such method is that of behavioral assays: that is, acceptance/rejection of non/ nestmates. Such assays, however, need to be designed according to the aggressiveness level of the species. In species that are highly territorial and aggressive, simple dyadic encounters between nestmates and nonnestmates will result in an unequivocal aggressive response, while in less aggressive species similar assays may result in mild responses, irrespective of whether the encountering ants are nestmates or alien (Buczkowski, Kumar, Suib, & Silverman, 2005).
In the latter case, assays that determine whether two nests merge into a single colony may provide more unequivocal results (Boulay, Katzav-Gozansky, Vander Meer, & Hefetz, 2003). To date, only a few studies have employed multiple methods to assess social and population structures (Ellis, Procter, Buckham-Bonnett, & Robinson, 2017). In the present study, we employed all three methods: genetic (microsatellite analysis), chemical (CHC nest profiles), and behavioral (colony insularity/openness tests), and demonstrate that none of the methods when applied individually is sufficient in itself for delineating colony boundaries.
The genus also expresses various reproductive strategies, from sexual to the asexual reproduction of new gynes (Darras, Leniaud, & Aron, 2014;Eyer, Leniaud, Darras, & Aron, 2013;Leniaud, Darras, Boulay, & Aron, 2012). Cataglyphis is distributed throughout the Old World deserts and arid habitats. A recently identified new species, C. israelensis, is distributed mainly in the north and north-eastern parts of Israel (Eyer, Seltzer, Reiner-Brodetzki, & Hefetz, 2017;Ionescu & Eyer, 2016), and is the subject of our study. The study was conducted at the distribution edge of this species along the northern coastline of Israel, presenting a similar ecological niche to that of other recently studied congener species (Leniaud et al., 2011;Saar, Leniaud, Aron, & Hefetz, 2014) that display various forms of sociality.

| Sample collection
Desert ants in the genus Cataglyphis are diurnal, highly thermophilic ants that are wide spread in arid habitats of the Old World. Their nest demography and architecture seem to be constrained by ecological factors. For example, in habitats along the Mediterranean Sea shores, nest depths are limited compared to the much deeper nests located in inland areas (Clémencet & Doums, 2007). Cataglyphis are central-place foragers, usually scavengers, but also contribute to pollination and seed dispersal (Boulay, Carro, Soriguer, & Cerdá, 2007;Cerda, Retana, & Cros, 1997). The ants forage individually and do not use any form of recruitment other than performing simple invitation behavior when a food source is discovered (Amor, Ortega, Cerdá, & Boulay, 2010). They are also well known for their high navigational ability, using both path integration and landmark guidance mechanisms (Wehner, 2008). Like many other ant species, they commonly have well-insulated nests and use CHCs as nestmate recognition cues (Hefetz, 2007;Soroker et al., 1994). They present highly diverse reproduction strategies, social structures, and population structures (Boulay et al., 2017;Lenoir et al., 2009  Only ten of these nests contained queens, always a single queen per nest. Because the Habonim population is sparse, and sample size was low for this population, it is included in only some of the analyses. Number of nests in each analysis is provided in Supporting Information Table S1.

| Genetic analysis
Since the study area hosts at least two sympatric Cataglyphis species in the bicolor group, we first confirmed that all the nests collected indeed belong to C. israelensis by genotyping Cytochrom B using one worker per nest (Eyer et al., 2017). CytB primers used: CB1 (Forward) TATGTACTACCATGAGGACAAATATC and CB2 (Reverse) ATTACACCTCCTAATTTATTAGGAAT, annealing temperature was 48°C (Simon et al., 1994). PCR products were purified with the spin-column PCR purification kit (TIANGEN Biotech), followed by sequencing with the ABI BigDye Terminator v.3.1 Cycle Sequencing Kit (Applied Biosystems, Foster City, CA, USA). Sequencing was performed on an ABI 3,730 Genetic Analyzer (Applied Biosystems). Base calling and sequence reconciliation were inspected by eye and performed using CodonCode Aligner (CodonCode Corporation, Dedham, MA, USA). Sequences were aligned using MUSCLE algorithm (Edgar, 2004) implemented in CodonCode Aligner, together with sequences previously used in Eyer 2017. All samples that were used in this project had the same CytB mitotype, similar to the C. israelensis mitotype from the main geographic distribution of the species.
Genetic (microsatellites) analysis was performed on 348 workers from the Atlit population (12*29 nests) and 72 workers (12*6 nests) from the Habonim population. DNA was extracted with 5% CHELEX (BIO-RAD) and amplified with eight microsatellite markers that had been designed previously for C. hispanica (Ch23, Ch08, Ch01, Ch10, Ch11 (Darras, et al., 2014)), C. cursor (Cc99, Cc54 (Pearcy, Clémencet, Chameron, Aron, & Doums, 2004), and C. niger (Cn02 (Saar et al., 2014) using Type-it PCR mix (QIAGEN). PCR products were sequenced by ABI3500 sequencer, and genotypes were then visualized with GeneMarker. Microsatellite data were checked for linkage disequlibrium, fixation index, and G-test using GENEPOP (Raymond & Rousset, 1995), heterozygosity was checked using both GENEPOP and FSTAT (Goudet, 1995), and F st was checked using FSTAT. To assess the number of queens, queens' genotypes, and the number of matings we used the program COLONY (Jones & Wang, 2010). The program STRUCTURE (Pritchard, Stephens, & Donnelly, 2000) was used to evaluate the number of colonies in the population. The number of K (number of possible colonies, from 1 to n size) was evaluated by the Evanno method using STRUCTURE Harvester (Earl, 2012). PCA of microsatellite data was performed as described by Ryan et al. (2017;using R version 3.3.3; library adegenet and ggplot2 for visualization), and relatedness was estimated using Coancestry (Wang, 2011; Queller GT method).

| Cuticular hydrocarbons analysis
For cuticular hydrocarbon analysis, we immersed individual heads in hexane (5-10 ants from each of 36 nests), containing 100 ng/ F I G U R E 1 (a) From left to right: study site locations at Atlit and Habonim (7 km apart), and nests sampled at the Atlit site: colors of the nests are according to COLONY matrilines (same colors indicate the same colony as in table ST2), colors of the nests are according to STRUCTUE analysis. (b) STRUCTURE harvester results indicating most probable K number. (c) STRUCTURE results indicating K = 6 or 22; individuals are grouped by nest origin, and each assumed colony is color coded. (d) PCA of microsatellite data with colony colors as indicated by STRUCTURE µl of C20 as internal standard. Initial analysis was conducted by gas chromatography/mass spectrometry (GC/MS), using a VF-5 ms capillary column, temperature-programed from 60 to 300°C (with 1 min initial hold) at a rate of 10°C per min, with a final hold of 15 min. Compound identification was achieved according to their fragmentation pattern, and the respective retention indices were compared to authentic compounds. Quantitative analyses were performed by flame ionization gas chromatography (GC/FID), using the above running conditions. Peak integration was performed using the program Galaxie Varian 1.9. Only hydrocarbons (a total of 24 compounds), which presumably play a major role in nestmate recognition (Lahav, Soroker, Hefetz, & Vander Meer, 1999), were used in the analyses. The relative amounts of each peak were calculated as a percentage of all the compounds in the analysis.
Hydrocarbon profiles of individuals were first plotted using principal component analysis (PCA; prcomp function of stats plotting was done by autoplot in ggplot2); however, this unsupervised analysis does not provide much information when the differences between samples are small. We then used k-means clustering in order to understand the structuring of the population, and assessed the most probable K using k-means function in stats, as well as mclust in mclust and nbclust in NbClust libraries in R. Plot was done by clusplot function in cluster library. Hydrocarbon profile specificity was also determined using multivariate statistics, that is, discriminant analysis using the stepwise forward mode, using IBM SPSS statistics 21. This module of the program works initially as a principle component analysis for reducing the number of variables to fit the number of cases and then performs the discriminant function.
We used the colonies suggested by the STRUCTURE program for discrimination.

| Behavioral assays
For the nest merging experiments, we composed seven queenright (QR) colonies comprising approximately 200 workers that were maintained in nesting boxes prior to the experiment. Using various combinations, we conducted merging assays between two colonies, which were chosen either according to the dissimilarity of CHC profiles or to nest metric distance (Table 1). We also performed merging assays between QR colonies and their neighboring QL nests, with each of the three QL nests being introduced to all of the QR nests, consecutively. Workers of each colony were marked with different colors using nail polish. At the onset of the merging experiment, the two nests were connected to a common foraging arena (size 20 × 20 cm). Sand and food were placed in the foraging arena, while water and cotton balls were placed in the enclosed nests. Scan observations were performed while registering interactions between nonnestmates, including aggression, trophallaxis, adult carrying, and antennation. Observations for the first four hours were carried out for 10 min at 10 min intervals, after which the nests were inspected every 12 hr for up to a week (Boulay et al., 2003). Merging tests between QR and QL nests were performed for the duration of only two hours for each test, with an interval of one day between each test involving the different QR nests.

| Association between genetic, chemical, and geographic distances
The pattern of nest distribution was assessed using the mean distance between nearest neighbors equation as described in Clark and Evans (1954). Pearson correlation was carried out between all types of distances (library SYNCSA, library matrix for data preparation, library ggplot2 and ggpmisc for visualization in R version 3.3.3). Mantel test was calculated by using the mantel.rtest function (library ada in R version 3.3.3) between all types of distances.
To test the hypothesis that the population is polydomous, we performed, using library dendextend, cluster analyses using each of the following distances: geographical (calculated from GPS locations), genetical (F st 's) and chemical (Mahalanobis distances). We then checked for association between the methods, in order to understand how they help to elucidate the population structure, by using cluster analysis (method: average), using the hclust algorithm in R version 3.3.3. Since each method resulted in a specific and slightly different dendogram, we assessed their validity by calculating the entanglement level between pairs of dendograms. Entanglement, a function used in the library dendextend, is measured by assigning values of consecutive numbers to the left dendogram's labels and then matching these numbers to those assigned to the right dendogram. Entanglement is the L-norm distance between these two vectors. In other words, we take the sum of the absolute difference (each one in the power of L; for example, sum(abs(x-y) L ) and divide it by the "worst case" entanglement level (e.g., when the right dendogram is the complete reverse of the left dendogram). L determines the penalty level (positive number, usually between 0-2).
L > 1 means that we give a big penalty for sharp angles; while L-> 0 means that whenever something is not a straight horizontal line, it receives a large penalty If L = 0.1, it means that we give higher preference to horizontal lines over nonhorizontal lines. We used L = 1.5 to allow for some flexibility; using lower L values results in higher entanglement value. A lower entanglement coefficient corresponds to a better alignment.

| Chemical analyses
Analyses of the hexane head extracts revealed a complex profile comprising 24 identified hydrocarbons (Figure 2a), accompanied by several branched fatty acids, alcohols, and esters (not shown in the chromatogram). Using the relative abundance of the different hydrocarbons, we first performed a principle component analysis that gave very little information, due to the small differences between samples (Supporting Information Figure S1). We then used the k-means approach to determine the level of structuring in the samples. The most probable number of K's was two according to all three analysis methods used, with the two components explaining only 51% of the variability (Figure 2b). Discriminant analysis was also performed while assigning the nests to colonies according to the results obtained by the STRUCTURE program. As depicted in Figure 2c, separation between colonies was significant (Wilks Lambada p < 0.05; function 1 explains 65.9% and function 2 16.1% of the variance, respectively).
There was, however, also a clear distinction between sites. Colonies five and six, which were very distinct from all other colonies, originated from the Habonim population; while those from the Atlit population were clumped together, but still distinct from each other.

| B EHAVI OR AL A SSAYS
Nest merging was carried out between QR nests of differing proximities (8-100 m; Table 1). In the assays involving nests that were far apart (90 and 100 m, respectively) workers showed only mild mutual aggression, and even entered each other's nests without any obvious reaction by the resident workers. However, on the third day the respective queens of both nests were executed, and a few days later also 10 to 20 of the workers (from each nest) were killed, indicating that the nests belonged to different colonies. In nests that were spatially nearer to each other (15 meters or less) we observed workers entering each other's nests and noted mutual high aggression between workers and from workers to queens (10-20 worker-worker attacks in 10 min); aggression was also noted in the common foraging arena. Because the queens were heavily attacked, we stopped the tests before they were injured.
In the assays involving QR (nests 100, 101, and 102) and neighboring QL nests (nests 103, 104, and 105) all the QL nests merged only with QR nest 102, which was located 2, 5, and 11 meters away from them, respectively.

| D ISCUSS I ON
The newly-described Cataglyphis israelensis is distributed along the central mountain ridge in Israel, but has also a small extension along the coast, at the foot of Mount Carmel, where it is sympatric with C. savignyi. In this study, we determined the social and population structure of C. israelensis in its coastal distribution. The results of nest excavations indicate that this species is monogynous. This was also supported by the genetic analysis of several nests where all workers' genotypes matched that of their mother queen. However, in some of the nests, workers' genotypes did not correspond with the assumed queen's genotypes. Two possible explanations can account for this genotype discrepancy. Some of the nest may have had multiple queens, which we may have missed during excavation. Alternatively, these workers may represent drifting workers that entered and were accepted in the wrong nest. The phenomenon of drifting workers is also known to occur in other Cataglyphis species (M. Knaden, personal communications), as well as other social insects (wasps, Sumner, Lucas, Barker, & Isaac, 2007;bumblebees, Lopez-Vaamonde, Koning, Brown, Jordan, & Bourke, 2004;and bees, Bordier, Pioz, Crauser, Conte, & Alaux, 2017). In a viscose population with very similar CHC profiles as we observed in our case, it could be elevated due to higher chance for individual-level mistakes. Drifting of workers may be also adaptive, being associated with mechanisms to increase within nest genetic diversity and a form of social parasitism (Nonacs, 2017;Zanette, Miller, Faria, Lopez-Vaamonde, & Bourke, 2014). On top of that, the queens are polyandrous, mating with five males on average. This concurs with findings for most of the Cataglyphis species studied to date (Boulay et al., 2017). Nest abundance in the studied populations was 0.00185 nests/m 2 , with each nest comprising 200 to 2,000 workers. To unravel the social and population structure of C. israelensis at the study site, we employed three different methods.
First, genetic microsatellite analysis (STRUCTURE) indicated a highly polydomous population, and possibly a supercolony. However, although the mean F st was rather low (0.15), it was not as low as that of a previously described supercolony of C. niger (0.005; (Leniaud et al., 2011), but more similar to that described for the monodomous populations of C. savignyi and C. livida (0.14 and 0.21, respectively (Timmermans et al., 2010;Leniaud et al., 2011) Similar results were obtained from the cuticular hydrocarbon analyses, with most of the nests, mainly from the Atlit population, demonstrating similar profiles, as expected from a unicolonial structure. Despite this profile similarity, the discriminant analysis was still able to discriminate several colonies. However, some of the nests that exhibited similar profiles were geographically remote from each other (and thus cannot be considered to share the same colony, unless they are part of a very large supercolony). CHC profiles are influenced both genetically and environmentally. Thus, in a small area like the one examined here, the possibility of isolation by distance (as suggested from the genetic data) superimposed by similarity in prey composition, may make differences between colonies hard to detect. Another variable to consider is that in our analysis we examined only the hydrocarbon constituents (following their demonstrated role in nest insularity in the closely related species C. niger (Lahav et al., 1999), but we cannot exclude the possibility that in this species polar compounds too are used for colonial identity.
The third method for delineating colony boundaries used in this study was to examine behavioral nest insularity/openness. Although, when given the opportunity, the presumed allocolonial ants entered each other's nests without demonstrating overt mutual aggression, they did not really merge into one large colony but, rather, invaded the opposite nest and killed its respective queen. Hence, these results refute any option of the existence of a supercolony structure, as may have otherwise been interpreted from the genetic or chemical analyses alone. The results of the experiments in which the queenless nests readily merged with a nearby queenright nest are consistent with a polydomous population structure.
To clarify some of the ambiguity in the results obtained by the different methods used in the study, we carried out a comparison among these methods. We first clustered the nests into assumed colonies according to either their physical location or microsatellite composition, or CHC profiles, and then determined whether these colony assignments were aligned. The underlying assumption for such a comparison is that the population structure is polydomous and that polydomous nests should be located in close proximity to each other. Neither the cluster based on microsatellites profiles nor that based on CHC profiles gave complete alignment with the cluster based on distances between nests. Furthermore, the finding that the microsatellite analysis did not match that of the CHC analysis indicates that CHC composition is not influenced by genetic composition alone. The comparisons further indicate that each method, when used individually for colony delineation, gives slightly different results. However, we cannot exclude the possibility that our markers may not have been sensitive enough to determine population structure of this small and somewhat isolated and viscous population.
In conclusion, the most probable population structure of C. israelensis along the coastline of Atlit is that of polydomy, with each colony consisting in one main nest, in which the single queen is located, and featuring several additional satellite nests, apparently queenless. This population structure seems to be very common in this genus and habitat, especially in the bicolor group (C. savignyi; Saar et al., 2014, C. drusus;unpublished data). The basis for such population structure could be due to difficulty of dispersal, properties of the substrate, and/or richness of the habitat. When isolation by distance occurs, new gynes hardly disperse and daughter colonies are located in high proximity to mother colonies. In soft substrate like sand, it is easier to dig new nests and, therefore, nest turnover may be large. Finally, a high availability of food increases nest mass and is perhaps conducive to colony spread.
Our study also emphasizes the importance of using a multimethod approach in population structure studies. A multi-method approach is often used for delimiting species' boundaries and we suggest its use also in delimiting colony boundaries. Most of the methods frequently used to determine a specific colony signature are affected by factors such as isolation by distance, food abundance and distribution, social structure, and mating systems. Thus, using a single method for assessment is insufficient and the importance of behavioral assays cannot be overemphasized.

ACK N OWLED G M ENTS
We would like to thank the reviewers for their considered com-