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

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.

INTRODUCTION

Understanding how animal societies are organized constitutes a major question in evolutionary biology (Wilson, 2012).Comparing social and population structures can shed light on the evolution of sociality as well as on other traits (Kennedy et al., 2017). In social Hymenoptera, although social structures and reproductive strategies can vary greatly among species (Boomsma, 2009; Crozier & Pamilo, 1996), 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, 1964, 1972 ; Hö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 Cataglyphis is considered a good model for social evolution (Boulay et al., 2017; Lenoir, Aron, Cerda, & Hefetz, 2009). It comprises over a hundred species (Agosti, 1990), displaying various social and population structures, ranging from the basic monogynous multicolonial population to polygynous supercolonies (Jowers et al., 2013; Leniaud, Heftez, Grumiau, & Aron, 2011; Timmermans, Grumiau, Hefetz, & Aron, 2010; Timmermans, Hefetz, Fournier, & Aron, 2008). 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.

METHODS

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).

Nests of C. israelensis were collected from two distinct populations on the coastal plain of northern Israel, in semi‐stabilized sand dunes, just south of Mount Carmel: from Atlit (N33 04.607 E35 06.554) comprising 37 nests in a plot of 100 × 200 m; and from Habonim (N32 38.920 E34 55.669; about 7 km south of Atlit; Figure 1a) comprising eight nests in a plot of 100x350 m. Nests were thoroughly excavated to include the queen and all workers, and brood. To ensure that the entire nest was collected, excavation was further extended by an extra 0.5–1 m in all directions after having observed the last gallery. For genetic and chemical analysis, all the nests in the experimental plot were excavated over four consecutive days.

Figure 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

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/µ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.

Table 1

List of nests used in the various merging experiments. Behavior occurrence is indicated by +or—(if the behavior did not occur). QR‐queenright nest, QL‐queenless nest

Nest A	Nest B	Merging	Between‐colony trophallaxis	Between‐colony aggression	Maximal distance between colonies (m)
		(+/−)	(+/−)	(+/−)
7 QR	21 QR	−	−	−	100
8 QR	26 QR	−	−	−	90
101 QR	100 QR	−	−	+	15
100 QR	102 QR	−	−	+	14
101 QR	102 QR	−	−	+	8
100 QR	103/104/105 QL	−	−	+	18,24, and 16
101 QR	103/104/105 QL	−	−	+	7,11, and 3
102 QR	103/104/105 QL	+	+	−	2,11, and 5

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)) 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.

RESULTS

Genetic analyses

None of the eight microsatellite markers used in this study showed indications of null alleles or linkage disequilibrium. The number of alleles ranged from 4 to 29. Mean heterozygosity was H o = 0.65 (range 0.44–0.69) H e = 0.107 (range 0.065–0.128). The fixation index Fit was slightly higher than zero (0.001 ± 0.03 SE jackknife) indicating random mating in the tested population. F st was very low at 0.15, while the within‐nest relatedness was 0.3 (±0.019 SE jackknife), which is considerably lower than assumed for a monogynous, monandrous nest (0.75). Nest excavations, in completion, revealed a single queen in nine nests. In 11 nests, COLONY genetic data analyses showed that all of the workers genotypes matched to a single queen, confirming unequivocally monogyny. In the other nests, some of the workers’ genotypes (average of 2.79 workers per nest; SE = 0.26) did not perfectly match the assigned queens’ genotype. These workers were observed both in QL (queenless) and in QR (queenright) nests; noteworthy, these QR nests had a single queen. Similar reconstructed matrilines are assumed to belong to the same colony (Supporting Information Table S2 and in Figure 1a). COLONY analysis also indicated polyandry, with an average of five males contributing to worker production (range 3–8; SE = 0.28).

Population analysis using the program STRUCTURE indicated that the population consisted in six or 22 colonies out of the 35 nests analyzed (Figure 1a–c). In Habonim the six nests sampled constituted two or three colonies, while in Atlit the 29 nests sampled constituted of four or 20 colonies. The latter value (20 colonies in Atlit) seems less probable in light of the nonuniform genotype distribution per nest, and the fact that the STRUCTURE program assigned different colonies for individuals that had been collected from the same nest. Thus, using this method alone is problematic in such cases. Similar results were produced using PCA (Figure 1d), which enables observation of the similarity of the genotypes between the groups.

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.

Figure 2

(a) Typical cuticular hydrocarbon profile of C. israelensis. The chromatogram shows only the 46 identified hydrocarbons (numbered peaks), of which 24 hydrocarbons were used in the analysis (marked in bold). (b) K‐means of CHC profiles K = 2, individuals are marked with nest number (first 3 digits) and individual ID. Circle indicates k‐1, triangle k‐2. The two components explain 51.02% of the variation. (c) Discriminant analysis based on CHC profiles while pooling nests to their assigned genetic colony (colored as depicted in Figure 1c), as indicated by STRUCTURE (k = 6)

BEHAVIORAL ASSAYS

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.

Figure 3a,b depicts the correlation between geographic and chemical or genetical distances, and Figure 3c the correlation between the chemical and genetical distances. There was no correlation between geographic and chemical distance (R 2 = 0.014; Pearson p value = 0.03; Figure 3a). In contrast, there was a positive correlation between the genetic F st's and the geographic distance between the nests (R 2 = 0.37; Pearson p value ≈ 0 Figure 3b), and a slight positive correlation between the chemical and genetic distances (R 2 = 0.19; Pearson p value ≈ 0 Figure 3c). Mantel test gave comparable results: there was low correlation between chemical and geographical distances (r = 0.124, p = 0.115), and positive correlation between genetic with geographic distances (r = 0.427, p≈0), as well as between chemical and genetic distances (r = 0.39, p≈0).

Figure 3

Regression analyses between cuticular hydrocarbon profiles (a) and F st values (b), respectively, and metric distances between nests. Red line indicates linear regression, gray indicates 95% confidence interval (R 2 = 0.014 and 0.37 respectively). (c) Correlation between CHC and genetic distances

Calculation of the nearest‐neighbor extrapolated‐nest‐density revealed an R value of 0.419 (R = mean distance of nearest neighbor observed/expected), indicating that the nests are aggregated more than expected under random distribution. We therefore checked for nest clusters using cluster analysis based on average geographic distances between nests (GPS locations). In Atlit there were eight main clusters of nests that were less than 20 meters apart from each other. Using the same nests, we performed a second cluster analysis based on average genetic distances (F st's). In this clustering, we obtained a slightly larger number of clusters (12 main clusters; Supporting Information Figure S2), with F st < 0.05. Aligning the two dendograms into a tanglegram resulted in a high level of entanglement (E = 0.63 in Atlit only), indicating a poor match between the dendograms. A third cluster analysis with the same nest population was performed using the average chemical Mahalanobis distances identified nine clusters (x < 15), but nest composition of each of these clusters differed from that based on geographic distances, as revealed by the high‐entanglement value (0.75). Finally, comparison of the dendograms based on genetic and chemical distances similarly resulted in a mismatch with an entanglement value of 0.64.

DISCUSSION

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/m2, 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)). Moreover, these F st values, along with our finding that they are positively correlated with distance between nests, can also be interpreted as a result of isolation by distance. We suggest that because this population is at the species’ distribution edge, as well as having fragmented distribution due to increasing land use for construction etc., or to different ecological niches, nests are somewhat isolated and gynes have difficulty dispersing. Alternatively, the genetic gradient with distance could be the consequence of their mating strategy, such as dispersal by budding, cloning of queens, or simply reduced flight ability of the gynes, as described for other Cataglyphis species. It follows that unraveling only the genetic structure of colonies provides a rather ambiguous estimate of population structure. There may be several explanations of why the STRUCTURE program collapsed nests into a unicolonial structure. Within‐nest relatedness was rather low (0.3) due to the high number of matings per queen, but the F st values still differed from zero, indicating a genetic separation between nests. Thus, colonies can be closely related but nonetheless do not form a unicolonial structure. This might also be due to a reasonably high level of drifting workers. Comparing nests with similar matrilines (using COLONY) seems to give a more accurate and biologically meaningful result, however, because the method of STRUCTURE is commonly used, the comparison is noteworthy.

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 multi‐method 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.

AUTHORS’ CONTRIBUTIONS

TRB and AH conceived the ideas and designed the methodology; TRB collected the data; TRB analyzed the genetic and behavioral data; TRB and AH analyzed the chemical data; TRB and AH led the writing of the manuscript. All authors contributed critically to the drafts and gave final approval for publication.

DATA ACCESSIBILITY

Sampling locations, CHC data and microsatellite genotypes are available at: Dryad https://doi.org/10.5061/dryad.6gb3977.

Click here for additional data file.

Click here for additional data file.

Click here for additional data file.

Click here for additional data file.