Reproductive variation. We first examined reproductive variation among individuals associated with collection from inside versus outside of the nest in the AZ samples and associated with seasonal changes (early versus late season collection) in the CA samples. From AZ, individuals collected from inside of the nest had significantly smaller ovary size measures than individuals collected outside of the nest (P < 0.001, Mann-Whitney U-test, Fig. 1). Individuals collected from inside of the nest compared to those collected outside of the nest did not significantly differ in the frequency of yellow bodies found on the ovaries (χ2 = 1.76, P = 0.18), though ovary size was significantly larger in individuals with yellow bodies, which indicate previous reproduction (P = 0.005, Mann-Whitney U-test, Fig. S1). In the CA population, individuals collected in the early season (April-July) had significantly larger ovary size than individuals collected in the late season (August-November) (P < 0.001, Mann-Whitney U-test, Fig. 2). In addition, significantly fewer individuals from the early season were found to have yellow bodies on their ovaries compared to late season individuals (χ2 = 4.95, P = 0.026). Individuals with yellow bodies were found to have smaller ovaries than individuals without yellow bodies (P = 0.02, Mann-Whitney U-test, Fig. S2), indicating possible post-reproductive individuals with small ovary size. Across all samples, we found that the ovary size index was significantly correlated with the Dufour’s gland width (R2 = 0.43, P < 0.001) and the Dufour’s gland length (R2 = 0.06, P = 0.035). Comparing datasets from the two populations, the individuals from California had a greater ovary size index (CA mean = 1.05, AZ mean = 0.54, df = 68.36, t = -4.78, P < 0.001) with significant differences in variance between populations (Levene’s test, P < 0.001). Intertegular distance was also significantly different between populations (CA mean = 6.81, AZ mean = 7.30, df = 75.25, t = 4.13, P < 0.001), though variances were not significantly different (Levene’s test, P = 0.97).
Chemical variation. We detected 42 total compounds across extraction types, after filtering as described in the methods. Twenty-nine compounds were found in the cuticular extracts and 32 compounds were found in the Dufour’s gland extracts. The cuticular and Dufour’s gland extracts showed clear separation in NMDS chemical space (Fig. 2), having both quantitative and qualitative differences, with 19 compounds shared between them. Both the cuticular and Dufour’s gland extracts consisted of alkenes and alkanes, with fatty acid esters also present in both, and one alcohol found on the cuticle (Table 1). The cuticular extracts contained multiple short chain alkenes and alkanes not found in the Dufour’s gland extracts, which in turn had several heavier fatty acid esters not found in the cuticular extracts. In general, the cuticular extracts had a greater percentage of the overall profile consisting of lighter molecular weight compounds when compared to the Dufour’s gland (Fig. 3). The Dufour’s gland extracts in X. sonorina are broadly comparable to Dufour’s gland extracts from two Xylocopa species characterized by Williams et al., 1983, consisting primarily of alkenes, alkanes, and esters. No other characterization of female Xylocopa cuticular chemistry has been conducted, but the detected cuticular compounds are also comparable to those seen in other bee species (Kather & Martin, 2015).
Considering only the AZ population, individuals sampled inside versus outside of the nests had significantly different cuticular profiles in PERMANOVA analysis (F1,29 = 3.60, P = 0.015) with no differences in dispersion (F1,29 = 0.25, P = 0.618). According to variable importance analysis, two probable alkenes (pentacosene, retention time: 24.9; octacosene, retention time 33.2) were most responsible for differences in cuticular extracts between individuals sampled from inside versus outside of nests (Fig. S3). Dufour’s gland extracts between individuals collected outside of the nest versus the inside of the nest were also significantly different (F1,16 = 2.70, P = 0.048) with no differences in dispersion between the groups (F1,16 = 0.67, P = 0.425). Pentacosane (retention time: 25.3) appeared to contribute the most to these differences (Fig. S4).
In the CA population, the cuticular extracts were not significantly different between individuals collected in the early vs late season (F1,39 = 0.67, P = 0.63), though there was a significant difference in their dispersions (F1,39 = 5.71, P = 0.021). In addition, Dufour’s gland extracts were not significantly different between early and late season sampled individuals (F1,23 = 2.01, P = 0.10) with no difference in dispersion between the two groups (F1,23 = 0.63, P = 0.43).
The chemical profiles between the two populations overlapped substantially in both the cuticular and Dufour’s gland extracts (Fig. 2), with no quantitative differences found between them. Given the two separate collection approaches for the CA and AZ dataset and the differences in ovary size/variance between datasets, we sought to minimize the effects of this on assessing qualitative population differences by comparing individuals from the AZ and CA populations that were likely to be at a similar point in their seasonal and life-history cycles. Consequently, to tentatively compare qualitative differences between populations, we compared AZ individuals that were collected outside of the nest to CA individuals collected in the early season (excluding three individuals from CA collected from the inside of a nest). The expectation is that this set of individuals from both CA and AZ should consist primarily of foraging individuals from the breeding season. The cuticular extracts between populations for this set of individuals (AZ, n = 18; CA, n = 20) were not significantly different in PERMANOVA analysis (F1,37 = 1.89, P = 0.098) with no significant difference in their dispersion (F1,37 = 1.57, P = 0.22). Similarly, for this set of individuals, the Dufour’s gland extracts (AZ, n = 11; CA, n = 12) were not significantly different between populations (F1,21 = 1.76, P = 0.12) with no significant differences in dispersion (F1,21 = 0.54, P = 0.47).
Correlations between chemical and reproductive variation. We examined correlations between all compounds and ovary size separately by population and extraction type, using FDR correction to account for the total number of individual correlation tests run (n = 122). Full results of correlation analysis, with Spearman’s rho values and corrected and uncorrected p-values can be found in Table S2-S3, though we provide an overview here.
Considering the cuticular extracts, in the AZ population, 17 of the 29 compounds were significantly correlated with ovary size, representing an average of 28% of the total relative abundance of compounds in the profile. In the CA population, 8 of 29 compounds were significantly correlated with ovary size, representing an average of 44% of the total relative abundance of compounds in the profile. Six compounds were significantly correlated with ovary size in both populations. In the CA samples but not AZ samples, ovary size was significantly correlated with heptacosane and nonacosane, compounds of roughly 20% and 11% relative abundances, respectively, across populations (Table 1). This is the primary reason why the AZ samples have more individual compounds sensitive to ovary size but show a lower overall percent relative abundance of ovary sensitive compounds.
Considering the Dufour’s gland extracts, in the AZ population, 17 of 32 compounds were significantly correlated with ovary size, representing approximately 39% of the total relative abundance of compounds in the profile. In the CA population, seven of the 32 compounds were significantly correlated with ovary size, representing approximately 11% of the total relative abundance of compounds in the profile. All seven of the compounds in the CA population that were correlated with ovary size were also correlated with ovary size in the AZ population.
Although strong correlations were found with ovary size in both extraction types and in both populations, only one compound, a putative isomer of pentacosene, was significantly correlated with ovary size across both the cuticular and Dufour’s gland extracts in both populations. However, all the compounds within each extraction type that showed significant correlation with ovary size in both populations also shared the same basic relationship (positive or negative).
Finally, we sought to better understand possible reasons for the differences in the degree of correlation between chemistry and ovary size across the two populations. While CA early season bees and AZ bees collected from outside of logs show similar chemical profiles and likely represent bees from similar behaviors/life stages, the CA late season bees and AZ bees from inside of logs are potentially very different in their life-history characteristics, which may influence the strength of correlations observed. In PERMANOVA analysis, the cuticular extracts from AZ bees collected inside of the nest (n = 13) compared to CA bees collected late season (n = 19) were significantly different (F1,31 = 2.86, P = 0.024), with no differences in dispersion (F1,31 = 0.576, P = 0.454). The Dufour’s glands of AZ bees collected inside of the nest (n = 7) compared to Dufour’s glands from CA bees collected late season (n = 11) were significantly different in PERMANOVA analysis (F1,16 = 3.22, P = 0.029), with no differences in dispersion (F1,16 = 0.056, P = 0.816).
Table 1. Relative percent abundances and standard error measurements for compounds in the cuticular and Dufour’s gland extractions for each population. Bolded cells with asterisks indicate a significant correlation with ovary size index in that extraction/population combination following FDR correction. “ND” indicates that the compound was not detected in that extraction/population combination. Compound order corresponds to the arrangement of compounds by retention time as shown in Fig. 4. For the “compound ID” column the letter “C” followed by a number indicates an alkane or alkene of that expected carbon length. Possible structural isomers of alkenes are indcated by “_1”, “ _2”, or “ _3” following the carbon length number in the compound ID column.