Mitochondrial and nuclear diversity of colonies of varying origins: contrasting patterns inferred from the intergenic tRNAleu-cox2 region and immune SNPs

Abstract In this study, we gathered sequence data from the tRNAleu-cox2 intergenic mitochondrial (mtDNA) region concurrently with single nucleotide polymorphism (SNP) data from 91 loci of nuclear DNA (ncDNA). The data was obtained from 156 colonies sampled in six apiaries from four countries. The full dataset was analysed and discussed for genetic patterns with a focus on cytonuclear diversity and admixture levels.

Introgression; SNPs; mitochondrial DNA; cytonuclear diversity; tRNA leu -cox2 haplotypes The honey bee, Apis mellifera L., natural variation is represented by over 24 subspecies belonging to four lineages: African (A), Middle eastern (O), Western European (M), and Eastern European (C) (Ruttner, 1988). The long-term sustainability of A. mellifera is largely dependent on locally adapted populations. However, during the last decades, globalization stimulated the introduction of commercial breeds in many regions of the A. mellifera native range which have hybridized with the local populations .
In this work, the cytonuclear background of honey bee colonies sampled in four countries was assessed using the tRNA leu -cox2 mitochondrial region and 91 ncDNA immune SNPs. A total of 156 samples were collected from single apiaries in Portugal, Spain, France, and Israel (Figure 1;  Supplementary Table S1). In France, the samples were from Avignon (N ¼ 24) and the island of Ouessant (N ¼ 12). In Portugal, they were from Braganc¸a (N ¼ 23) and the island of São Miguel,Azores (N ¼ 23). Spain and Israel were represented by samples collected from an apiary in Marchamalo (N ¼ 34) and in Beit Dagan (N ¼ 40), respectively. Each sample was composed of a pool of 10 workers collected from a single colony. Total DNA was extracted from pools of 20 forelegs using Ron's Tissue DNA Mini Kit (Bioron).
The tRNA leu -cox2 region was PCR-amplified with the primers E2 and H2 using the reaction conditions of Garnery et al. (1993). PCR products were Sangersequenced by STABVIDA (Portugal). The alignment of the sequences was performed with MEGA 6.06 (Tamura et al., 2013). Sequences of two samples were discarded from the mtDNA analysis because they failed quality standards.
Twenty haplotypes (27 variants) were identified in the 154 samples, of which 11 (13 variants) belonged to lineage A, seven (10 variants) to lineage M, and two (four variants) to lineage C (Supplementary Table S1). One haplotype (M81; accession number MZ767556), as well as eight variants, are reported for the first time (Supplementary Table S2). This novel variation was detected in Portugal and Spain, further supporting the view of Iberia as a hotspot of  Table  S1). In contrast, in Spain, 82% of the colonies belonged to lineage M of which 53% carried M4, the most common haplotype in northeastern Iberia (Ch avez-Galarza et al., 2017). These results are consistent with earlier maternal surveys in Iberia, reporting a cline formed by haplotypes of M-lineage in the northeastern part of Spain and A-lineage in Portugal and southwestern part of Spain and virtually no Clineage introgression in the native Apis mellifera iberiensis (Ch avez- Galarza et al., 2017;Garnery et al., 1998;Miguel et al., 2007). In São Miguel, despite previously reported high levels of C-lineage introgression, a single A14 haplotype (the novel A14c variant) was detected in the sampled colonies and this is the most common and widespread African haplotype in the Azores (Ferreira et al., 2020).
In France, all colonies from the Avignon apiary were of C-lineage ancestry (91% C2 and 9% C1), contrasting with the apiary of Ouessant where all colonies carried M-lineage mitochondria (58% M17e and 42% M4). These results are consistent with reports of high C-lineage introgression into the native Apis mellifera mellifera in many regions of mainland France and with Ouessant acting as a refuge for this subspecies (Miguel et al., 2007;Pinto et al., 2014;Requier et al., 2019). As for Avignon, in the apiary of Beit Dagan, Israel, all colonies were of C-lineage ancestry, with haplotype C2 also predominating (95%) over C1 (5%). While Apis mellifera syriaca (lineage O) is the native subspecies of Israel (Ruttner, 1988), our data is consistent with the claim that it has largely been replaced by commercial stocks of C-lineage ancestry (Soroker et al., 2018). However, according to anecdotal reports, the C-lineage A. m. ligustica has been favoured by Israeli beekeepers (Soroker et al., 2018) and the finding of mostly C2 in Beit Dagan suggests that most colonies in this apiary likely derived from Apis mellifera carnica or Apis mellifera caucasia imports, similar to Avignon. This assertion stems from the fact that while C2 is the predominant haplotype in the distributional ranges of those two subspecies, A. m. ligustica is characterized by harbouring haplotype C1 (Franck et al., 2000;Kandemir et al., 2006;Muñoz et al., 2009;€ Ozdil et al., 2009;Pinto et al., 2014). We caution, however, against any over-interpretation of the data presented herein as they were collected from single apiaries and might therefore suffer from biases linked to management practices such as colony splitting or queen rearing.
The 156 pooled samples were genotyped in the MassARRAYV R MALDI-TOF platform using a recently developed assay containing 91 immune SNPs (Henriques et al., 2021). The membership proportions (Q-values) were inferred from these SNPs with ADMIXTURE (Alexander et al., 2009) using the same settings and reference populations of Henriques et al. (2021). GENEALEX 6.50 (Peakall & Smouse, 2006) was used to estimate the mean number of haplotypes, the effective number of haplotypes, and unbiased diversity whereas allelic richness was calculated with ADZE v1 (Szpiech et al., 2008).
In the ADMIXTURE analysis, the most likely number of clusters was K ¼ 3 (cross-validation test), coinciding with lineages A, M and C detected with mtDNA in the four countries. At K ¼ 4, it is possible to distinguish not only lineages A and C (the purple and orange cluster, respectively), but also the two M-lineages subspecies A. m. mellifera (green cluster) and A. m. iberiensis (blue cluster) represented in Figure 1 (see Supplementary Table S1 and Figure 5 in Henriques et al., 2021, for details).
The membership proportions inferred from the SNP data in ADMIXTURE are concordant with the mtDNA data in what concerns C-lineage introgression. Colonies from Portugal, Spain, and Ouessant were virtually free of C-lineage cytonuclear introgression. In contrast, colonies from Avignon and Beit Dagan, which only harboured C-lineage haplotypes, exhibited membership proportions into the C-lineage cluster >32% (Figure 1, Supplementary Table S1). However, these colonies also displayed higher levels of genetic diversity for ncDNA but not for mtDNA, which contrasts with the other mainland colonies from Iberia showing the opposite trend (Table 1).
The finding of elevated nuclear diversity in admixed populations relative to parental populations is consistent with predictions from population genetics theory and has been confirmed by empirical data in honey bee studies (Harpur et al., 2012). Of note is that this finding is also a by-product of the biased nature of the SNP assay employed here (Henriques et al., 2021), which was tailored to include immunerelated SNPs segregating between the two most divergent lineages: M and C. This means that both alleles will co-occur in admixed populations whereas only one alternative allele is expected to be found in parental populations.
Despite the limited geographical sampling undertaken in each country, this study further corroborates the (i) genetic integrity of A. m. mellifera in Ouessant and A. m. iberiensis, (ii) high levels of mitochondrial diversity in Iberia, supported by the detection of a novel haplotype and many variants, and (iii) presence of C-lineage introgression in the native range of M-lineage (mainland France) and O-lineage (Israel) subspecies. Moreover, this study detected contrasting levels of diversity in the cytonuclear compartments. The higher nuclear diversity in the colonies from mainland France and Israel can be explained by admixture. While admixture may lead to increased diversity (Harpur et al., 2012), as measured by the classical parameters (Table 1), it may also lead to disruption of gene complexes that have been fine-tuned by selection acting on local variation across micro-evolutionary time scale . This type of diversity may be important for the populations to face the challenges of a rapidly changing world. Table 1. Diversity estimates for mtDNA (left column) and SNPs (right column) for each sampled apiary from Portugal (PT), Spain (SP), France (FR) and Israel (IL). Na -mean number of alleles; Ne -number of effective alleles; uhunbiased diversity; Ar -allelic richness.