Introduction
Pest management has played a critical role in achieving the current level of the global food supply, and its importance will continue to increase due to the increase in the human population with the Sahelian countries being particularly affected. In parallel, concern about preserving the environment and the health of living organisms is also increasing. In this context, biological control (BC) is often considered as one of the most economical and environmentally sustainable means of managing agricultural pests and crop diseases (Cock et al., Reference Cock, van Lenteren, Brodeur, Barratt, Bigler, Bolckmans, Cônsoli, Haas, Mason and Parra2010; Oleke et al., Reference Oleke, Manyong and Mignouna2013; van Lenteren et al., Reference Van Lenteren, Bolckmans, Köhl, Ravensberg and Urbaneja2018). BC of insect pests relies on the use of natural enemies (mainly parasitoids and predators, and microbial pathogens) which can be purchased from commercial suppliers or produced locally, or/and which occur naturally in the agroecosystem concerned. In augmentative BC (ABC), natural enemies are mass-reared before being released for supplementary pest control in field crops either for immediate or several generation control. ABC has been applied with some success for more than 100 years in several cropping systems and geographic areas (Gurr & Wratten, Reference Gurr and Wratten2000; Cock et al., Reference Cock, van Lenteren, Brodeur, Barratt, Bigler, Bolckmans, Cônsoli, Haas, Mason and Parra2010; Calvo et al., Reference Calvo, Bolckmans and Belda2012; Ba et al., Reference Ba, Baoua, Kaboré, Laouali, Oumarou, Dabire-Binso and Sanon2014 ; Wang et al., Reference Wang, He, Zhang, Lu and Babendreier2014; Kaboré et al., Reference Kaboré, Ba, Ba, Dabire-Binso and Sanon2017; van Lenteren et al., Reference Van Lenteren, Bolckmans, Köhl, Ravensberg and Urbaneja2018). Today, ABC is most often included in Integrated Pest Management (IPM) programmes that also incorporate new technologies and practices, in particular genetic tools and approaches (van Lenteren et al., Reference Van Lenteren, Bolckmans, Köhl, Ravensberg and Urbaneja2018). Genetic markers, including mitochondrial DNA and microsatellite loci, are particularly valuable for solving taxonomic problems, assessing advantageous biological traits such as host specificity and dispersal ability in natural enemies, detecting among- or within-population structure and diversity, etc. (Behura, Reference Behura2006; Gariepy et al., Reference Gariepy, Kuhlmann, Gillott and Erlandson2007; Rauth et al., Reference Rauth, Linz, Gerber and Hufbauer2011; Gaskin et al., Reference Gaskin, Bon, Cock, Cristofaro, De Biase, De Clerck-Floate, Ellison, Hinz, Hufbauer and Julien2011; Vorsino et al., Reference Vorsino, Wieczorek, Wright and Messing2012; Moffat & Smith, Reference Moffat and Smith2015). Concerning ABC, there are several important issues which require the use of molecular tools including species identification, detection of signatures of bottlenecks, loss of genetic diversity and increased inbreeding (the last two likely due to the small numbers of BC agents released), monitoring the released BC agents, and gene flow with natural populations in the agroecosystems targeted. All these parameters can impair BC efficacy and need to be better understood to improve our knowledge of what actually happens in the agroecosystems.
Pearl millet Pennisetum glaucum (L.) is one of the most important staple cereal foods in semi-arid and arid areas of Asia and Africa, and particularly in the African Sahel where it is the only cereal crop that can produce grain yields under the prevailing harsh conditions (i.e. the poorest soils, the hottest climate, low and erratic rainfalls). Despite its adaptation to extreme climatic conditions, pearl millet suffers from many biotic constraints including insect pests (Nwanze & Harris, Reference Nwanze and Harris1992). Among these, the millet head miner (MHM), Heliocheilus albipunctella (de Joannis) (Lepidoptera, Noctuidae), is a major insect pest of millet in the African Sahel. It is a univoltine species which diapauses during the dry season (i.e. for at least 9 months) and attacks millet during the rainy season (Gahukar et al., Reference Gahukar, Guèvremont, Bhatnagar, Doumbia, Ndoye and Pierrard1986; Nwanze & Sivakumar, Reference Nwanze and Sivakumar1990; Sow et al., Reference Sow, Brévault, Delvare, Haran, Benoit, Cœur d'Acier, Galan, Thiaw, Soti and Sembène2018). Management of this pest in Sahelian countries has involved the use of insecticides, host-plant resistance, cultural management practices and more recently the use of indigenous BC agents (Payne et al., Reference Payne, Tapsoba, Baoua, Ba, N'Diaye and Dabiré-Binso2011; Ba et al., Reference Ba, Baoua, Kaboré, Laouali, Oumarou, Dabire-Binso and Sanon2014; Kaboré et al., Reference Kaboré, Ba, Ba, Dabire-Binso and Sanon2017; Baoua et al., Reference Baoua, Ba, Amadou, Kaboré and Dabire-Binso2018). Habrobracon hebetor Say (Hymenoptera, Braconidae), a well-known gregarious ectoparasitoid of the larvae of a wide range of economically important moths that mainly infest stored grains as well as field crops worldwide is commonly used. This cosmopolitan parasitic wasp has been used with some success for many years to regulate destructive moth pests (Payne et al., Reference Payne, Tapsoba, Baoua, Ba, N'Diaye and Dabiré-Binso2011; Adarkwah et al., Reference Adarkwah, Ulrichs, Schaarschmidt, Badii, Addai, Obeng-Ofori and Schöller2014; Ba et al., Reference Ba, Baoua, Kaboré, Laouali, Oumarou, Dabire-Binso and Sanon2014; Ghimire & Philipps, Reference Ghimire and Philipps2014). The first releases of H. hebetor against MHM were attempted in Senegal in 1985, then in Niger in the early 2000s (Kaboré et al., Reference Kaboré, Ba, Ba, Dabire-Binso and Sanon2017). In recent years, the apparent efficacy of H. hebetor augmentative releases against MHM populations conducted at the beginning of the rainy season has been reported (Payne et al., Reference Payne, Tapsoba, Baoua, Ba, N'Diaye and Dabiré-Binso2011; Ba et al., Reference Ba, Baoua, Kaboré, Laouali, Oumarou, Dabire-Binso and Sanon2014; Baoua et al., Reference Baoua, Ba, Amadou, Kaboré and Dabire-Binso2018). As a result, its use has expanded but many knowledge gaps remain particularly concerning the genetic diversity and structure of this species that occurs naturally in open-field cropping systems and grain stores as well as being mass-produced and released in the field. For example, it has been commonly assumed that even a restricted parasitoid population can survive the 9 month diapause period of its host MHM. H. hebetor can also enter reproductive diapause but only under laboratory, and certain temperature and photoperiod conditions which do not occur in Sahelian countries such as Niger (Chen et al., Reference Chen, Zhang, Zhu and Throne2012). Consequently, as a parasitoid of many lepidopteran species, is it likely that H. hebetor shifts to other moth hosts and/or locations and recovers rapidly in the rainy season. The fact is, what actually happens to this widely used BC agent before (i.e. during mass rearing) and after its release among all the H. hebetor which occur in pearl millet fields is not clearly understood. To this end and to advance our knowledge on H. hebetor population genetics, a large sampling of H. hebetor mainly emerging from MHM larvae parasitized in pearl millet fields in Niger, one of the major pearl millet producers in Africa, was collected and compared to some laboratory strains. Based on two mtDNA genes and microsatellite markers (Garba et al., Reference Garba, Loiseau, Benoit and Gauthier2016), and using population genetics and phylogeny approaches, we investigated: (1) the genetic diversity of some laboratory H. hebetor strains used for annual field releases, (2) whether it would be possible to monitor them in the field, (3) the patterns of genetic diversity and structure of H. hebetor from the field, (4) the possible occurrence of recent bottlenecks and inbreeding, which are more likely to occur in the dry season when its host MHM diapauses and more broadly in the absence of many lepidopteran hosts, (5) the level of gene flow between the samples and the distinct genetic (sub)populations, and at various spatial scales. Our findings are discussed in light of the species’ biology and their relevance for BC.
Material and methods
Samples
A total of 618 H. hebetor individuals from 34 samples comprising seven laboratory samples and tests, and 27 samples collected in millet fields from 12 different locations in Niger (fig. 1) and Iran was studied (table 1). Thus, one sample corresponds of all individuals collected on the same day from several panicles of millet plants growing in a restricted site, or from the laboratory colony or tests. The samples named L-NgF, L-Ng9, L-Ng12, L-Ng14 and L-Ng19 consisted of specimens collected between August 2014 and September 2016 from the colony of H. hebetor maintained at the DGPV (Department of plant protection) in Niamey since 1989. This laboratory colony was first established from parasitized MHM larvae collected in the field in Sadoré (45 km from Niamey). It was subsequently maintained on Corcyra cephalonica 3rd and 4th instar larvae as described by Bal (Reference Bal2002). Both insects are cultured in the laboratory at room temperature (30 °C ± 4 °C) with 30% relative humidity (RH). To maintain the colony, H. hebetor from parasitized MHM larvae collected in the field or caught using insect nets were sometimes added, as was the case in 2008 using specimens from the Dosso region, and in 2010 and 2015 using specimens from Tondibia and Louga regions, respectively (fig. 2). The samples named C-Ng17 and C-Ng18 were the progeny of preliminary reciprocal crossing-tests between H. hebetor from the DGPV colony and some field individuals from Bazaga (GPS coordinates 13.79999N, 5.09999E_ about 325 kms east from Niamey) as follows: 2 ♂ field x 2 ♀ DGPV for C-Ng17, 2 ♀ field x 2 ♂ DGPV for C-Ng18. The other H. hebetor samples were from MHM larvae collected from either a restricted area including 16 sites along a maximal 1400 meter-transect (i.e. Louga), or from sites located between 2 and 500 kms away (fig. 1). In all the areas, millet was the most important cereal grown but was often associated with cowpea (Vigna unguiculata, Fabaceae). In the part of the Louga region sampled, there were temporary ponds 1–2 km away from the central collection point that enabled vegetable cultivation, and where the natural vegetation included the two predominant tree species Faidherbia albida (Fabaceae) and Piliostigma reticulatum (Fabaceae). As far as we know, some of the sites in which H. hebetor came from were previously subjected to some DGPV-parasitoid colony releases during BC campaigns or experiments: in Dantiadou, one release was made by the DGPV and a private producer in October 2014; in Dosso, two successive releases by the DGPV were made in the autumn 2014 and 2015; and others in the localities of Louga and Tondikouarey in September 2014, and in the locality of Ayorou in September 2015 (fig. 2). In the present study, all samples collected in the Louga region (except F-Ng36 and F-Ng37) resulted from a preliminary release experiment of parasitoids conducted by the DGPV in September 2016 in a millet field located at the central point (GPS coordinates 13.64946N, 2.09554E) on the 5 September 2016. Before the parasitoids were released at the central point, parasitized MHM larvae were collected in the vicinity and kept in the laboratory until H. hebetor emergence (i.e. F-Ng20, control sample). After this control collection, a laboratory wasp strain was released and dispersed to parasitize MHM larvae in millet fields. One week later, on the 12 September, parasitized MHM larvae were collected from nine millet panicles per site (F-Ng21 to F-Ng35) along cardinal directions (N, W, S), and at intervals of 200 m away the release point along a maximal 0 to 1400 m transect. Because of a severe drought and the vicinity of a transport axis, no panicles and/or MHM larvae were found along the eastern and southern axes. A total of 10–15 parasitized MHM were collected at each site, transferred to the DGPV laboratory, and examined for parasitoid emergence. The wasps obtained were clearly female biased. The wasps originating from MHM collected per panicle were pooled and placed alive in a tube filled with 90% ethanol and stored at −20 °C for genetic analyses. From one to 30 wasps emerging from different MHM larvae collected from different panicles to limit sibling effects per site were then analyzed. Most of these specimens were also morphologically checked by the authors and when required by a Braconid expert (G. Delvare, CIRAD, France).
Niger maps (from OpenStreetMap) with the sites where parasitized millet head miner (MHM) larvae from which the H. hebetor samples emerged were sampled from September 2013 to September 2016 (fig. 1a). Ng-Lou represents all 12 samples from Louga collected in September 2016 (fig. 1b). The sites are from pearl millet fields. The sample code for each site is as listed in table 1.
Figure representing what is know about the contributions of field H. habrobracon in the DGPV laboratory colony (from 1987 to 2015) and the releases of H. hebetor from the DGPV colony (autumn 2014 and 2015) in some areas where the sampling investigated was conducted (from September 2013 to September 2016).
Description of the samples and number of Habrobracon hebetor analyzed (N). Samples were collected from millet crop fields in Niger and Iran, and from DGPV laboratory rearings or obtained from a cross between specimens sampled in the field and those reared in the laboratory.
Louga (T0 or T7- cardinal axis distance from the release point O) T0: the 5 September, 2016; T7: the 12 September, 2016; N: North; W: West; S: South.
DGPV: Department of Plant Protection, Niamey, Niger.
In front of the sample code, the following letters F for field, L for laboratory and ,C for crossing tests. Country names and sample codes are in bold.
DNA protocols
Genomic DNA was extracted from each adult H. hebetor for further amplification of mitochondrial and nuclear DNA markers using the DNeasy Blood & Tissue Kit (Qiagen) according to the manufacturer's instructions.
Two mitochondrial fragment genes, COI and 16S rDNA, were amplified and then sequenced in one to five specimens (male and/or female) per site and laboratory strain (table 1) using the barcoding primers LCOI490 [GGTCAACAAATCATAAAGATATTGG] and HCOI2198 [TAAACTTCAGGGTGACCAAAAAATCA] for COI (Folmer et al., Reference Folmer, Black, Hoeh, Lutz and Vrijenhoek1994), and the primers 16SWa [CGTCGATTTGAACTCAAATC] and 16SWb [CACCTGTTTATCAAAAACAT] for 16S rDNA (Dowton et al., Reference Dowton, Austin and Antolin1998). Standard 25 µl polymerase chain reactions (PCRs) were performed. For COI, 1.2 µl of one DNA template, 2.5 µl of Qiagen amplification buffer 10X (with 15 Mm MgCl2), l μl of 2.5 Mm dNTPs, 1.75 µl of each 10 µM primer, 1 unit of Taq polymerase (Qiagen), and 16.3 µl of ultrapure water were used. PCRs were performed as follows: an initial denaturing step for 5 min at 95 °C, 40 cycles comprising a 3 min. denaturing step at 95 °C, a 1 min. annealing step at 50 °C, a 1 min. elongation step at 72 °C, and a final 10 min elongation step at 72 °C. For 16S rDNA, 2.5 µl of one DNA template, 2.5 µl of Qiagen amplification buffer 10X (with 15 Mm MgCl2), l.8 µl of 2.5 Mm dNTPs, 2 µl of each 10 µM primer, 1 unit of Taq polymerase (Qiagen), and 13.7 µl of ultrapure water were used. PCRs were performed as follows: an initial denaturing step for 4 min at 95 °C, 37 or 40 cycles comprising a 30 s denaturing step at 95 °C, a 1 min. annealing step at 55 °C, a 40 s elongation step at 72 °C, and a final 7 min elongation step at 72 °C. Eighty-two amplified COI- and 78 16S rDNA fragments were sequenced (Eurofins MWG, France) and obtained in both directions. The resulting consensus sequences were aligned using Muscle software (Li et al., Reference Li, Cowley, Uludag, Gur, McWilliam, Squizzato, Park, Buso and Lopez2015; http://www.ebi.ac.uk/Tools/msa/muscle/) and manually checked before being deposited in GenBank (Accession no. MH733498-576 and MH733582-584 for COI; MH721433-435 and MH721438-513 for 16S rDNA), and analyzed.
A total of 614 H. hebetor individuals originating from fields, laboratory strains or crossing tests was first screened using 24 microsatellite loci distributed in three PCR multiplex sets (table S1). Standard 10 µl PCRs were performed as described by Garba et al. (Reference Garba, Loiseau, Benoit and Gauthier2016). Diluted PCR products were run on an ABI Prism 3130XL automated sequencer (Applied Biosystems, Montpellier, France) using the GeneScan-500 LIZ ™ size standard. Allele sizes were scored using GeneMapper ™ 4.0 software (Applied Biosystems) and confirmed manually.
Genetic analyses
Hardy–Weinberg and linkage testing, marker polymorphism analysis
GenePop 4.2 on the Web (Raymond & Rousset, Reference Raymond and Rousset1995, http://genepop.curtin.edu.au/) was used to compute the observed number of alleles (N A), the allele size range and mean observed (H O) and unbiased expected heterozygosity (H E) of all samples containing more than five specimens. The same software was also used to test linkage disequilibrium (LD) between each pair of microsatellite loci, deviation from Hardy–Weinberg equilibrium (HWE) and the inbreeding coefficient (F IS) of each locus and sample. As multiple tests were conducted, sequential Bonferroni correction of the P-values was performed. As potential factors of deviation from HWE, null allele frequency and/or the presence of scoring errors were estimated (Micro-Checker program, van Oosterhout et al., Reference Van Oosterhout, Hutchinson, Wills and Shipley2004; FreeNA package, Chapuis & Estoup, Reference Chapuis and Estoup2007). Standard indices of genetic variability for the COI- and 16S rDNA fragments obtained [number of haplotypes (N H), polymorphic sites, haplotype (h) and nucleotide (π) diversities and their standard deviations, nucleotide difference] and the selective neutrality of the observed nucleotide polymorphisms estimated using Tajima's D parameter (Tajima, Reference Tajima1989) were evaluated using the DnaSP v6.11.01 package (Rozas et al., Reference Rozas, Sanchez-Delbarrio, Messeguer and Rozas2003).
Phylogenetic inferences
Phylogenetic COI- and 16S rDNA-based trees were built using character-based maximum likelihood (ML) methods (Guindon et al., Reference Guindon, Dufayard, Lefort, Anisimova, Hordijk and Gascuel2010; PhyML3.0 http://www.atgc-montpellier.fr/phyml/). The Smart Model Selection tool (Lefort et al., Reference Lefort, Longueville and Gascuel2017) in the PhyML Web server was used to determine the most suitable model of DNA substitution for each dataset. Based on the Akaike information criterion (AIC), the GTR + I and GTR + G + I models were selected as the statistically most appropriate models for the 16Sr DNA and COI sequence datasets, respectively. Branch supports were assessed by 1000 bootstrap resamplings. Phylogenetic trees were constructed with 94 or 93 sequences from the present study: 82 COI- and 78 16S rDNA H. hebetor sequences, 12 COI- and 15 16S rDNA sequences from various Braconids [Bracon sp1. and sp2., B. brevicornis, B. lissogaster, B. nigricans, B. cephi], some more H. hebetor retrieved from GenBank [KX371823, KJ627789, KY484509, AF003483] and the BOLD database [NZHYM803, GBMI73841], and rooted.
Population genetic structure
A multilocus Bayesian clustering analysis was performed using a Markov Chain Monte Carlo (MCMC) algorithm and Structure software v.2.3 (Pritchard et al., Reference Pritchard, Stephens and Donnelly2000) to delineate genetic populations. The most likely value of K (number of genetic populations), the quantification Q (individual assignment) of how likely each individual is to belong to each genetic population were estimated using the admixture model with correlated allele frequencies, and using (or not) information on their prior location (LocPrior = 1) (Hubisz et al., Reference Hubisz, Falush, Stephens and Pritchard2009). For each run, a burn-in period of 80,000 iterations followed by 80,000 MCMC iterations was applied. Log-likelihood estimates [lnP(D)] were calculated for K = 1–15. Twenty replicates were performed to test the consistency of the results. The most probable number of genetic populations (K) was estimated as the lowest value of K which best captured the structure observed (Structure v.2.3) and the higher delta K value estimated according to the method of Evanno et al. (Reference Evanno, Regnaut and Goudet2005) (Structure Harvester program, http://taylor0.biology.ucla.edu/structureHarvester/, Earl & vonHoldt, Reference Earl and vonHoldt2012). The dataset was then partitioned into each of the defined genetic populations to be reanalyzed following the same procedure as above. The values of r, which parameterizes the amount of information vehicled by the location, were also examined. Values of r > 1 indicate that there is either no population structure or that the structure is independent of the location. As a complementary analytical approach, principal component analysis (PCoA) was performed on the whole sampling and then separately on each genetic population to summarize the genetic variability of the microsatellite dataset and to identify possible genetic (sub)populations using GenAlex 6.5 (Peakall & Smouse, Reference Peakall and Smouse2012).
Levels of population genetic diversity
For each defined genetic (sub)population, the previously described genetic variability indices [(N A), (H O), (H E) for the microsatellites, (N H), (h) and (π) and their standard errors for the mtDNA markers] were estimated using GenePop on the Web and DnaSp, respectively.
As demographic bottlenecks can influence allele richness and heterozygosity, tests were run to identify their signatures in the microsatellite data for each of the 31 samples with N ≥ 5 using Bottleneck v.1.20.2 (Piry et al., Reference Piry, Luikart and Cornuet1999). This method assumes that recent bottlenecks cause a shift away from an L-shaped distribution of allelic frequencies to one with fewer alleles in low-frequency categories. So, we used the descriptor (‘mode-shift’ indicator) of the allele frequency distribution and the non-parametric Wilcoxon sign rank test under the two-phase mutation model (TPM) (variance 30.00, probability 70%, 1000 simulations) to assess deviation from expected excess heterozygosity (signature of a bottleneck). Probability values determined using a one-tailed Wilcoxon test for heterozygote excess <0.05 and a shifted-mode of allele frequency distribution were considered as indicative of bottlenecks.
Levels and factors of population differentiation
The levels of genetic differentiation among the 31 samples with N ≥ 5 and the (sub) populations delineated were assessed from each microsatellite dataset harbouring (1) no null alleles using the excluding null alleles (ENA) method (F ST-ENA) following bootstrapping over loci (1000 replications) using FreeNA (Chapuis & Estoup, Reference Chapuis and Estoup2007), and (2) null alleles (F ST) using GenePop 4.2. The overall significance of genotypic differentiation was estimated using Fisher's exact tests implemented in GenePop 4.2 on the Web.
The correlation between geographic distance and genetic distance was estimated by isolation-by-distance (IbD) analysis performed using a regression of F ST/(1 − F ST) values against the log of the geographical distance (kilometres) between sampling sites. The significance of the correlation between the two data matrices was assessed using a Mantel test with 1000 permutations, as implemented in Isolde in GenePop 4.2. on the Web. This approach was conducted at the local (i.e. Louga) and regional scales including all field samples with N ≥ 5.
To estimate the significance of genetic variation explained by the following factors: the collection year for the field samples (i.e. from 2013 to 2016) in Niger, the location of the field samples (i.e. Louga, Dantiadou, Tondikouarey, Bangoula, Tondibia, Lourgou, Ayorou, Ouallam, Dosso, Mali or Nigeria borders), the type of sample (lab vs. field samples), some analyses of molecular variance (AMOVA) were performed using GenAlex 6.5 (Peakall & Smouse, Reference Peakall and Smouse2012). The significance of each analysis was tested using 999 random permutations.
Results
Hardy–Weinberg and linkage testing, marker characteristics
A total of 614 H. hebetor specimens sampled across 34 samples were successfully genotyped at 24 polymorphic loci. The number of alleles per locus ranged from two to 24 (overall loci, mean N A = 8.91). Heterozygosities revealed a wide range of values per locus (0.018 < H O < 0.726; 0.012 < H E < 0.889) and high mean values over all 24 loci (H O = 0.516 + 0.205; H E = 0.643 + 0.206) (table S1). After Bonferroni correction, 21 loci out of 24 conformed to HWE and no evidence for LD remained, suggesting they segregated independently. The mean frequency of null alleles for each locus (table S1) and sample (table 2) was low (f < 0.10) and no evidence for significant scoring errors was revealed. As a result, the 21 microsatellite loci that conformed to HWE were used for the population genetics study. For all the samples with N ≥ 5, the number of alleles per locus ranged from 1.71 to 5, and mean heterozygosities (H O = 0.564 ± 0.085; H E = 0.533 ± 0.091) did not globally differ except for F-Ng15 and F-Ir12. But most of the samples individually revealed heterozygosity excess and low or even negative F IS (−0.442 ≤ F IS ≤ 0.133) values both statistically significant in 7/25 field samples, 1/2 crossing samples and 1/4 laboratory samples (table 2).
Genetic characteristics of the 34 samples over the 21 polymorphic microsatellite loci selected for the study (no LD, HW Equilibrium, f null−allele < 0.10).
*Samples with N < 5.
N, number of specimens genotyped; N A, mean number of alleles per locus per sample; H O, observed heterozygosity; H E, unbiased expected heterozygosity; SE, standard error.
f null−allele, null allele frequency (Chapuis & Estoup, Reference Chapuis and Estoup2007) significant f null−allele (f > 0.20).
F IS, inbreeding coefficient (Weir & Cockerham, Reference Weir and Cockerham1984).
Departure from Hardy–Weinberg equilibrium after Bonferroni correction are in bold (P < 0.05) with a test for heterozygote excess.
Twelve 641 bp-COI and 32 506 bp-16S rDNA fragments were obtained. Both genes revealed two very common haplotypes (f ≥ 15.8%) and many rare ones. No haplotype was shared between the H. hebetor from Iran and Niger. Specimens from laboratory samples shared haplotypes with those from the fields (table S2). The genetic diversity (h COI = 0.616 ± 0.057, π COI = 0.005 ± 0.002; h 16S rDNA = 0.584 ± 0.066, π16S rDNA = 0.004 ± 0.001) and the polymorphism (6.86% and 10.47% for COI- and 16S rDNA with <7% of informative sites, respectively) were shown to be moderate. This polymorphism did not result from non-random processes (Tajima's D COI = −1.976, P COI < 0.05; Tajima's D 16S rDNA = −2.413, P 16S rDNA < 0.01). It was mainly characterized by transitions among the H. hebetor haplotypes in Niger. The haplotype(s) deriving from Iran specimens exhibited a low level of nucleotide variability among themselves, in some cases, none, but a high level of both transitions and transversions with the homologous haplotypes deriving from Niger specimens (tables S2).
Phylogenetic inferences
Both COI- and 16S rDNA phylogenetic reconstructions (only the COI-based reconstruction is shown here, fig. 3) showed two well-supported H. hebetor lineages (bootstrap values >90%), differing in 5.7–6.3% of nucleotides in the 641 bp-COI fragment and 1.9–4% in the 506 bp-16S rDNA fragment. One lineage consisted of the haplotype(s) deriving from the Iran specimens. The other included all the haplotypes from the Niger specimens, displaying <1% of nucleotide variability for COI and between 0.3 and 2.7% for 16S rDNA (table S3). ML analyses of both mitochondrial datasets derived from H. hebetor from Niger failed to produce well-resolved reconstructions (most bootstrap values < 60%). Unexpectedly, all B. brevicornis haplotypes clustered within the H. hebetor lineage instead of being positioned as outgroups with the others braconids, and low nucleotide variability (1–1.3% for COI, 0.3–2.3% for 16S rDNA) was estimated between the two putative species (table S3, fig. 3).
Maximum likelihood reconstruction (Model GTR + G + I, γ = 0.411) showing phylogenetic clusterings among the 641-bp mtCOI sequences obtained in this study from H. hebetor specimens and other braconids (Bracon sp1, Bracon sp2, B. brevicornis, B. lissogaster, B. nigricans, B. cephi) as outgroups, homologous sequences from H. hebetor retrieved from Genbank and BOLD databases. The tree with the highest log likelihood (-Ln = 2124,53456) is shown. Bootstrap values over 60 are shown at the corresponding node. Abbreviations for H. hebetor [Hap for Haplotype 641 bp-COI, sample code acronym followed by the number of the specimen, or sample code followed by Hh, the gender (male M/female F) and the number of the specimen for the NgF sample] and for the other braconids are Brevi for B. brevicornis, Lisso for B. lissogaster, Cephi for B. cephi are the same as in table 1 or table S2.
[Hap-1]: F-Ng1-24, F-Ng1-4, F-Ng7-5, L-Ng12-10, F-Ng3-23, F-Ng6-3, F-Ng5-18, F-Ng7-13, F-Ng8-15, F-Ng5-2, F-Ng8-16, F-Ng4-13, F-Ng8-3, F-Ng6-11, F-Ng5-5, F-Ng3-2, F-Ng13-1, L-Ng14-12, L-Ng14-24, C-Ng17-1, C-Ng17-13, L-Ng19-8, F-Ng24-6, F-Ng24-9, F-Ng36b-22, F-Ng36-2, F-Ng37-1, F-Ng37b-1, F-Ng34-12, F-Ng34b-1, F-Ng34b-10, F-Ng21-5, F-Ng21-8, F-Ng23-11, F-Ng23-3, F-Ng25-6, F-Ng25-21, F-Ng28-1, F-Ng28-2, F-Ng28-4, F-Ng31-5, F-Ng31-13, F-Ng32-1, F-Ng33-16, F-Ng33-11, F-Ng35-13, F-Ng26-2, F-Ng26-9, F-Ng29-1.
[Hap-2]: F-Ng8-7.
[Hap-3]: F-Ng30-13.
[Hap-4]: F-Ng30-8.
[Hap-5]: F-Ng20-8, F-Ng20b-6.
[Hap-6]: F-Ng29-3.
[Hap-7]: F-Ng11-3, F-Ng11-18, F-Ng11-21, F-Ng34-5, F-Ng32-2.
[Hap-8]: F-Ng3-10.
[Hap-9]: L-Ng-HhF1, L-Ng-HhF2, L-Ng-HhM3, L-Ng-HhM4, L-Ng-HhM5, F-Ng1-23, F-Ng1-10, L-Ng9-21, L-Ng9-8, L-Ng9-24, F-Ng4-10, L-Ng12-15, F-Ng35-19.
[Hap-10]: F-Ng4-7, F-Ng13-6, F-Ng15-4, F-Ng15-5.
[Hap-11]: L-Ng9-9.
[Hap-12]: F-Ir1-2, F-Ir1-20, F-Ir1-12.
Population structure
Based on PCoA in which the first two axes explained respectively, 25.8 and 23.6% of the variance in the experimental data (data not shown), and Bayesian analyses (with or without LocPrior = 1), the H. hebetor from Iran were clearly separated from those from Niger and a complex structure was revealed within the Niger genetic population (K1) (fig. 4. fig S1). Using Evanno's method, the most likely number of subpopulations was seven (table S4), with most of the individuals being assigned (with Q ≥ 0.80) to one of the seven subpopulations shown in different colours and named K 1A to K 1G (fig. 4). Most individuals from Louga were more weakly assigned to a particular subpopulation (Q ≤ 0.70) with many individuals being genetically admixed (i.e. with K 1B, K 1C, K 1D, K 1E nuclear background) and shared among at least four subpopulations (K 1C, K 1D, K 1E, K 1F) (fig. 4). K 1A consisted of all individuals from L-NgF and L-Ng12, and some more from L-Ng9 and F-Ng1 individuals (0.60 < Q < 1). K 1B consisted of all individuals from L-Ng14, C-Ng17, C-Ng18, F-Ng13 and F-Ng15, and a few from F-Ng28. The H. hebetor included in K 1C were from L-Ng9, from most of the field samples collected at the regional and Louga scales in Niger but the probability of assignment to this subpopulation ranged from weak to high (0.35 < Q < 1). Admixture (major Q ≤ 0.50) was significant in many Louga samples (F-Ng24 to F-Ng29, F-Ng31 and F-Ng36-37) meaning some specimens were not assigned to a specific subpopulation. K 1D consisted of all L-Ng19 individuals released at Louga on the 5 September 2016 and 76% of the F-Ng21 individuals from MHM larvae collected in the same site but 1 week later. K 1E consisted of individuals from various locations (i.e. Tondikouarey, Tondibia and Louga). K 1F consisted of some individuals from two nearby sites in Louga. K 1 G only included the H. hebetor F-Ng11 from Dosso (table 3, fig. 4). The r values estimated (0.7 < r whole sampling < 1.6; mean r Niger = 0.468) indicated that the clear and complex structure shown between Iran and Niger and then in Niger, respectively, was not significantly based on location information.
Inference of population genetic structure in H. hebetor based on 21 polymorphic microsatellite loci. Results of successive Bayesian multiloci clustering analysis and the Evanno's method run on the samples from Iran and Niger, from the laboratory colony and cross-tests, and the field with the delineation of two main genetic populations (K1 in Niger vs. K2 in Iran) and seven subpopulations in Niger (K 1A-K 1G).
Each vertical bar represents one of the 614 specimens genotyped. The length of each bar is proportional to the inferred ancestry values (Q) in each genetic (sub)population for each specimen. The sample code for each sample is as listed in table 1. * for Ng36-Ng37.
Description and genetic diversity of the seven genetic subpopulations (K 1A to K 1G) defined in Niger, with the additional non-assigned specimens, using Structure software and Evanno's method.
The sample code for each sample is as listed in table 1 (i.e. L for Laboratory, F for field and C for crossing tests). The index 1 to 6 added to L, C and F indicates the chronological positioning of collection dates relative to each other (1: September 2013, 2: August 2014, 3: September 2014, 4: October 2014, 5: February 2015, 6: August 2015, 7: September 2015, 8: December 2015, 9: September 2016).
N. ass., not assigned specimens: when all Q ≤ 0.33.
N, number of specimens clustered; (p), proportion of specimens concerned per sample; N A, mean number of alleles per locus; H O, observed heterozygosity; H E, unbiased expected heterozygosity; SE, standard error.
Departure from Hardy–Weinberg equilibrium after Bonferroni correction are in bold (P < 0.05) with a test for heterozygote excess.
The sample code for each sample is as listed in table 1.
Patterns of genetic diversity
Overall, the highest levels of genetic diversity were obtained for the field samples vs. the laboratory ones for both markers [(microsatellites: H Ofield = 0.589 ± 0.119; H Otests = 0.534 ± 0.028; H Olab = 0.517 ± 0.044; mean N Afield = 3.77 ± 0.86; mean N Atests = 3.07 ± 1.24; mean N Alab = 3.16 ± 0.64), (COI mitotypes: three for laboratory specimens vs. 12 for the field ones)]. Among the field samples, significantly lower genetic diversity and 17/34 private alleles (i.e. alleles that occur at any frequency but in a single sample) were observed for F-Ir12 (table 2, table S2). In Niger, K 1C, which was also the most represented subpopulation (i.e. approximately 50% of the H. hebetor genotyped and many admixed individuals, most sites) displayed the highest level of nuclear (table 2) and mitotype diversity [9/12 COI haplotypes, 14/32 16S rDNA haplotypes, tables S2]. K 1F and K 1 G exhibited a marked excess of heterozygotes as previously reported for F-Ng32, F-Ng33 and F-Ng11 (table 2), and the lowest number of alleles observed (table 3).
Wilcoxon tests and the ‘mode-shift’ indicator implemented in the Bottleneck software indicated recent demographic bottlenecks for all the laboratory strains, C-Ng17 and some field samples [i.e. Tondikouarey (F-Ng3), Tondibia (F-Ng5), Mali and Nigeria borders (F-Ng13, F-Ng15), Louga (F-Ng20, F-Ng21, F-Ng24, F-Ng30)].
Patterns of genetic differentiation
F ST−ENA values (from FreeNA) and F ST values (from Genepop) were shown to be similar, which is consistent with the low frequency of null alleles found (table 2), so only the F ST−ENA values are presented here. A wide range of pairwise F ST−ENA values was observed (0.025 ≤ F ST−ENA ≤ 0.628), the highest values being with the H. hebetor from Iran and the lowest being between most of the samples from Louga (mean pairwise F ST−ENA (Ng23−Ng35) = 0.074 ± 0.031). This low level of genetic differentiation between 12 samples from Louga (F-Ng23 to F-Ng35), led us to pool them under the name ‘Lou’ even though Bayesian analyses separated ‘Lou’ H. hebetor into three different subpopulations (191 in K 1C vs. 12 in K 1E, 20 in K 1F). In contrast, the two central collection point samples from Louga, F-Ng20 (before the H. hebetor release, in K 1E) and F-Ng21 (after the H. hebetor release, six specimens in K 1C and 19 in K 1D), were shown to be moderately differentiated from the others (0.140 < F ST−ENA (Ng20−21) vs (Ng23−Ng35) < 0.226) and even more strongly differentiated among themselves (F ST−ENA Ng20 vs 21 = 0.304) (table 4, table S5; fig. 4). At the Niger scale, pairwise F ST−ENA values ranged from 0.055 to 0.507, the highest values were always observed with F-Ng15 (mean F ST−ENA = 0.373) and F-Ng11 (mean F ST−ENA = 0.299), which was apart in K 1 G (tables 3 and 4; fig. 4).
Pairwise estimates of F ST−ENA among the 31 samples from the laboratory colony, parasitized millet head miner larvae collected in millet fields and crossing tests over the 21 polymorphic loci. Specimens from F-Ng23 to F-Ng35 (Lou) sampled in the Louga region were pooled.
*Samples with < 5 specimens studied (N < 5). Sample names and types are in bold.
The laboratory samples showed levels of genetic differentiation between themselves that increased with the collection date from the laboratory colony: the pairwise F ST−ENA values between L-NgF, L-Ng9 and L-Ng12 (sampled from August 2014 to February 2015) were 0.087, with L-Ng14 (collected in December 2015) they ranged from 0.133 to 0.155 and then with L-Ng19 (collected in September 2016) they ranged from 0.191 to 0.291. Some were less differentiated from some field samples than from some laboratory strains (e.g. F ST−ENA (L−Ng14 vs F−Ng13) = 0.016, F ST−ENA (L−Ng14 vs F−Ng5) = 0.098 vs F ST−ENA (L−Ng14 vs L−Ng19) = 0.191). Regarding L- Ng19, all the estimated pairwise F ST−ENA values were greater than 0.19, with both reared and field samples. C-Ng17 and C-Ng18 which corresponded to reciprocal crossings between DGPV specimens and specimens collected in the field at Bazaga in September 2015, revealed low genetic differentiation from L-Ng14, moderate one from F-Ng1, F-Ng5, F-Ng13 and Ng-Lou samples and high genetic differentiation from all the others (table 4). The laboratory specimens thus clustered in four subpopulations (K 1A, K 1B, K 1C, and K 1D) (fig. 4) also including field specimens. Most comparisons between the F ST−ENA values were significant (P < 0.05). Regarding pairwise F ST−ENA values between each subpopulation defined in Niger (0.036 < F ST−ENA < 0.433), subpopulations K 1A, K 1D, K 1F, and K 1 G differed significantly from one another and from K 1B, K 1C, K 1E, K 1H (78.6% of the specimens studied) which, however, displayed low levels of genetic differentiation (F ST−ENA < 0.1) among them (table S6).
Factors of genetic differentiation
At local and regional scales, there was no significant evidence for isolation by distance (R 2 ≤ 0.09) between the field samples. AMOVA showed that genetic variation among the 24 field samples (i.e. with N > 5) or the 11 locations, or among the 24 field samples vs. the four laboratory ones (i.e. with N > 5) was not significantly explained by the collection year (3%), the location (5%), the type of sample (3%), or the samples (10–14%). The great majority of genetic variation was found within individuals (85–86%) which were consequently very distinct.
Discussion
This study is the first attempt to understand the patterns and the main driver(s) of population genetics in H. hebetor a biological agent of lepidopteran moths and a major pest of pearl millet fields in the Sahel. Consistent with previous work (Antolin et al., Reference Antolin, Ode, Heimpel, O'Hara and Strand2003), our findings highlight the importance of the mating system in the genetic diversity and structure of H. hebetor population. They also advance our understanding of its reported efficacy against MHM in pearl millet fields.
Main pattern of genetic diversity
All parameters estimated based on our mitochondrial and nuclear marker datasets showed moderate levels of genetic diversity across the H. hebetor sampling conducted in Niger and Iran but with a significantly lower level for the Iran sample, which also exhibited many private alleles and haplotypes. However, as only one sample was investigated, no conclusion can be drawn until other samples from this country are studied.
In Niger, the samples of H. hebetor emerging from MHM moths collected in millet fields displayed a higher level of genetic diversity than those sampled from the laboratory colony, which were actually collected at different times (i.e from August 2014 to September 2016). All laboratory samples unlike field ones were seen to have experienced a recent bottleneck. These findings are consistent with BC agents bred in mass-rearings before being released into agroecosystems. Indeed, mass-rearing generally begins with specimens collected in the field. The history of the laboratory colony of H. hebetor was partially known, but when a founder sample is collected in the field, the genetic variability of this base sample is assumed to be more significant, as demonstrated by our results. During rearing, the more adapted specimens to the rearing conditions, for example to their alternative host, will be maintained. If the effective number of breeding individuals decreases drastically at this stage, the laboratory colony might face a bottleneck and an increased inbreeding which will reduce the genetic variability of the reared strain before being restored by field specimens. Taken together, these phenomena could result in strains that are genetically less diverse, and genetically different from field populations. In the present study, contrary to what was reported in a previous work (Heimpel et al., Reference Heimpel, Antolin, Franqui and Strand1997) and as observed in other BC agents (eg. Streito et al., Reference Streito, Clouet, Hamdi and Gauthier2017), the artificial pressure exerted over many generations on the laboratory strain, did not allow the fixation of typical laboratory mitotype(s). Indeed they were also shared with field specimens indicating that the laboratory rearing had been ‘contaminated’ with field specimens either intentionally (to re-establish the laboratory colony) or involuntarily (specimens from lab experiments or from the surroundings), or that both specimen types moved and mated as suggested by the reciprocal cross-tests conducted in the DGPV. Consequently, monitoring of field-released H. hebetor conducted in the framework of BC programs would not be possible using mitotype markers.
General outbreeding and heterozygote excess for ‘field’ samples vs. laboratory strains: possible mechanisms involved
The F-statistics showed that F is estimates of 21 samples out of the 31 tested were negative, even highly negative (−0.442 < F is < 0.133, overall mean values F is = −0.077) indicating that there was significant heterozygote excess and more outbreeding than would be expected randomly relative to HWE. Particularly concerned are the field samples from Niger. Heterozygote excess in insect populations is not as common as heterozygote deficiency, and has never been examined in H. hebetor species. What specific characteristics and/or demographic events could explain this finding? Overdominant selection favouring heterozygotes, negative assortative mating, and bottlenecks are generally proposed as explanations for heterozygote excess in natural populations (Cornuet & Luikart, Reference Cornuet and Luikart1996; Greenspoon & M'Gonigle, Reference Greenspoon and M'Gonigle2014). Concerning overdominant selection favouring heterozygotes, the problem always consists of having convincing examples of gene loci at which the heterozygote was selectively superior. In H. hebetor no evidence for such selection has been reported in the literature. Thus, could negative assortative mating and/or bottleneck events explain the heterozygote excess situations observed in our study? Negative assortative mating, is nonrandom mating in which mating between individuals of the same type (e.g. same morph, same brood, etc.) occurs less often than expected. This nonrandom mating has mainly been demonstrated to contribute to the avoidance of inbreeding, the maintenance of polymorphism and to reducing the production of less fit offspring (Tien et al., Reference Tien, Massourakis, Sabelis and Egas2011; Rolán-Alvarez et al., Reference Rolán-Alvarez, Saura, Diz, Rivas, Alvarez, Cortés, de Coo, Estévez and Iglesias2012; Holman et al., Reference Holman, van Zweden, Linksvayer and d'Ettore2013; Jiang et al., Reference Jiang, Bolnick and Kirkpatrick2013; Follett et al., Reference Follett, Calvert and Golden2014; Greenspoon & M'Gonigle, Reference Greenspoon and M'Gonigle2014). The fact is that the mating behaviour and life history attributes of H. hebetor facilitate outbreeding. Wasp larvae pupate in close proximity to the consumed host but after emergence females exhibit a pre-mating refractory period during which both sexes disperse (Antolin et al., Reference Antolin, Ode, Heimpel, O'Hara and Strand2003; van Wilgenburg et al., Reference Van Wilgenburg, Driessen and Beukeboom2006). As H. hebetor is a gregarious species (i.e. more than one offspring can develop on/in a single host) the emergent wasps may be related and dispersal before mating consequently avoids inbreeding. Females are also able to distinguish males from their own natal brood thus avoiding related partners, and mate more than once and with different males (Antolin & Strand, Reference Antolin and Strand1992; Antolin et al., Reference Antolin, Ode, Heimpel, O'Hara and Strand2003). Negative assortative mating favouring outbreeding therefore occurs in this species. What is the particular advantage of outbreeding, and inherent heterozygote excess in H. hebetor? Like in many Hymenoptera species, sex determination is influenced by a genetic process called complementary sex determination (CSD) which, in H. hebetor, is based on a single sex-determining locus (sl-CSD). Under CSD, heterozygous diploids at the sex-determining locus are female, while homozygous haploid and diploid individuals at this locus develop into males. These diploid males generally comprise a genetic load within populations of Hymenoptera. In some cases, these males have limited viability or/and fail to pass their genetic material to subsequent generations because they are unable to mate properly or because they are sterile. Some diploid males have been shown to produce viable sperm, but the sperm is diploid rather than haploid and results in sterile triploid progeny (Cowan & Stahlhut, Reference Cowan and Stahlhut2004). To conclude, lack of fitness and low performance characterize these males. Because of the relationships among inbreeding, diploid male production, and reproductive failure by these diploid males, natural selection should favour outcrossing mechanisms in species like H. hebetor. Even though a female-biased sex ratio should be favoured to reduce the genetic load of the CSD both female and male-biased sex ratios and diploid males have been reported in natural populations (Antolin & Strand, Reference Antolin and Strand1992; Heimpel & de Boer, Reference Heimpel and de Boer2008; Antolin et al., Reference Antolin, Ode, Heimpel, O'Hara and Strand2003; Manishkumar et al., Reference Manishkumar, Dhirubhai and Piyushbhai2013; Ba et al., Reference Ba, Baoua, Kaboré, Laouali, Oumarou, Dabire-Binso and Sanon2014; Garba personal communication). Then, could bottleneck events also explain the heterozygote excess revealed in some of the field samples? During a bottleneck, as previously said, allelic diversity tends to decline more than heterozygosity because bottlenecks tend to eliminate many low-frequency alleles. The result is that there is an excess of heterozygosity in selectively neutral loci compared to what one would expect from a population at HWE. The heterozygosity excess persists only a few generations until a new equilibrium is established, consequently only recent historical bottlenecks can be detected. In H. hebetor, it has been assumed that the millet field parasitoid population would not survive the long off-season when its host, MHM, is in diapause, and consequently undergoes a serious reduction in its effective size and hence bottlenecks. Our analyses showed that only a minority of the field samples had faced significant recent bottlenecks but they did not show less genetic diversity than the other field samples. Interestingly, Tajima's negative values suggest population size expansion for most samples, which, in turn, no longer experienced bottlenecks. In species like H. hebetor whose generation time is <2 weeks under the climatic conditions of countries like Niger (Magro & Parra, Reference Magro and Parra2001; Eliopoulos & Stathas, Reference Eliopoulos and Stathas2008; Ba et al., Reference Ba, Baoua, Kaboré, Laouali, Oumarou, Dabire-Binso and Sanon2014) and can parasitize a wide range of lepidopteran hosts, in crop fields and in grain stores, recovery of a sufficient effective population size and genetic diversity may be rapid, at least in the field. Only one laboratory strain (L-Ng19), sampled in September 2016 from those studied showed a negative F is value similar to those observed in the field samples concerned by heterozygote excess and bottleneck. As discussed below, L-Ng19 specimens were sampled from laboratory colony only a few generations after some specimens from the Louga surroundings had been added to maintain the laboratory colony. This could have led it to differentiate from the previously sampled lab H. hebetor [F ST(L−Ng19 vs L−Ng−F−9−12−14) > 0.19] and to experience more outbreeding than those laboratory samples which had not been recently refreshed with field H. hebetor. Taken together, our results suggest that the mating system rather than bottlenecks is probably the main driver of genetic diversity in H. hebetor.
Complex genetic structure and gene flow
Insights from the two molecular markers were not equally significant, which may be explained by the features of the two kinds of markers; microsatellite loci, with their high mutation rate, can increase the overall resolution needed to detect more recent population genetics events, and in comparison, the slower mutation rate of the mtDNA markers may provide a more historical than recent picture of gene flow. Indeed, phylogenetic analyses based on COI- and 16srDNA datasets resulted in poorly resolved trees, with only two strongly supported divergent lineages distinguishing specimens from Iran from those from Niger. The lack of a phylogenetic signal within Niger H. hebetor clades may result from the low level of intraspecific divergence observed (e.g. <3% in mtCOI). As suggested by Mardulyn & Whitfield (Reference Mardulyn and Whitfield1999), the absence of a phylogenetic signal in parasitic braconid wasps could result from rapid diversification of the species examined caused by the great diversity of lepidopteran hosts. However, this hypothesis and obtaining a better phylogenetic signal could be tested by notably increasing the number of genes studied. Based on microsatellite loci, Iran and Niger H. hebetor were again clearly separated into two strongly differentiated populations (F st−ENA(Niger vs Iran) > 0.44). The extent to which all genetic parameters differed led us to wonder whether the two lineages were distinct morphocryptic entities of H. hebetor, possibly geographically associated as reported by Chomphukhiao & co-authors (Reference Chomphukhiao, Takano, Takasu and Uraichuen2018) or even two closely related species although both clustered with H. hebetor sequences. This question remains unanswered until additional samples from Iran can be analyzed.
Within the Niger H. hetetor lineage, microsatellite loci revealed a complex structure with seven subgenetic populations (K 1A to K 1G), most of which were a mixture of H. hebetor coming from distinct and even distant field locations and sites, laboratory strains, and different collection years. Neither the collection site and year nor the origin of the sample (laboratory vs field), but rather the individuals themselves significantly explained the observed genetic variation. No link could also be found between genetic and geographic distances for any field sample, indicating that human-mediated dispersal is significantly involved, which does not exclude the role of natural dispersal and other factors. This complex genetic structure was supported by a wide range of F ST−ENA values, expressing very diverse levels of genetic differentiation and gene flow both at local and regional scales in the millet agroecosystem, between as well as among laboratory strains and field samples, among and within the seven subpopulations.
Exactly what relevant insights does this genetic structure yield into the discriminatory power of the microsatellite loci, sample history, gene flows and H. hebetor dispersal? The subpopulation K 1A comprised the following laboratory strains: 100% of those sampled in August 2014 and February 2015, 61% of those sampled in October 2014, none of those sampled in December 2015 and September 2016, indicating that the colony has probably evolved over generations because of voluntary or/and involuntary ‘genetic contamination’, since some of them clustered with a remote (i.e. Dantiadou located about 70 km from Niamey) field sample. In addition, as far as we know, in late summer 2015, some H. hebetor from the vicinity of Louga were added to the DGPV colony, and the H. hebetor population(s) from this region were shown to consist of admixed and genetically distinct specimens therefore representing sources of new genetic diversity for the lab colony. Similarly, all the H. hebetor sampled in December 2015 among those reared in the laboratory (i.e. L-Ng14) clustered with the reciprocal crosses conducted between laboratory and field specimens in September 2015 and thus temporally closely related to L-Ng14, but also with all the H. hebetor sampled in two very distant field sites (i.e about 480 km from one another, about 350 km from the laboratory in Niamey). Here, it appears that human transport is the most likely explanation for the large-scale dispersal of this small parasitic wasp, since these sites are located far apart but host genetically similar specimens. The same conclusion can be drawn with the most represented H. hebetor subpopulation (K 1C) that encompassed most of the field specimens from Louga and elsewhere in Western Niger. However, human-mediated dispersal does not exclude the natural H. hebetor dispersal. Interestingly, all the specimens released at the central point in Louga (L-Ng19) clustered with 76% of those emerging from parasitized MHM larvae collected at the same site but 1 week later (F-Ng21) which suggests that these wasps might be the offspring of the released specimens a week earlier (L-Ng19). If so, in 1 week, the released H. hebetor did not disperse far from the point of release, <200 m, as no H. hebetor from this subpopulation was found farther away. Previous studies reported that H. hebetor dispersed 3–5 km from the release point in 5 weeks in Burkina Faso and Niger millet fields (Baoua et al., Reference Baoua, Ba, Amadou, Kaboré and Dabire-Binso2018) and 1 km in 6 days in a tomato field in Russia (Adashkevich et al., Reference Adashkevich, Saidova and Takanaev1986), significantly farther than we found in the current study. The evaluation of parasitoid dispersal capacity is an important parameter in the development of ABC programmes (Heimpel & Asplen, Reference Heimpel and Asplen2011). Further experiments should be conducted to evaluate the ability of this released-field parasitoid species to disperse and persist in the field and in grain stores during the rainy as well as in the dry season when hosts are scarce. Even though various markers and methods have been reported to mark and monitor insects (Hagler & Jackson, Reference Hagler and Jackson2001; Lavandero et al., Reference Lavandero, Wratten, Hagler and Jervis2004; Pomari-Fernandes et al., Reference Pomari-Fernandes, de Freitas Bueno, de Dortoli and Favetti2018), the microsatellite markers presently used might be also valuable for monitoring a ‘pure’ (i.e. not contaminated) laboratory strain released in field.
Interestingly, many admixture specimens were identified in Louga, making it difficult or even impossible to assign some of them to any subpopulation, thereby indicating high levels of gene flow. One possible reason for this particular observation may be that Louga is a crop and vegetable production area with shrubs and trees that could harbour alternative host insects and provide floral resources such as nectar for the H. hebetor population before or during the period when MHM larvae become available in millet fields. Additional vegetation also creates a suitable climatic environment for the parasitoids in the hot dry season. For example, H. hebetor wasps have already been observed in Faidherbia albida, a tree species commonly found in the Louga region. Consequently, connectivity between the ressources in this area was probably preserved in space as well as in time, at least during the rainy season, thereby favouring gene flow and diversity between the H. hebetor present. In accordance with this statement, 70% of the specimens included in the subpopulation that displayed the highest allelic richness and heterozygosity came from Louga. Then, when suitable hosts become scarce in the field they move to stored-grain product hosts. The F ST_ENA values between the Louga samples (F ST−ENA < 0.11) as well as the absence of isolation by distance reinforce these hypotheses by demonstrating a low level of differentiation presumably resulting from the homogenizing effect of the short-distance dispersal of the parasitoid between all the suitable hosts thereby driving a strong gene flow among the samples. On the other hand, the presence of some genetic differentiation among samples taken in the near vicinity (e.g. the release point vs. other samples at Louga), suggests that local management actions or environmental factors are unlikely to be countered by the natural dispersal of insects.
Surprisingly, F-Ng20 which represented the central point sample of Louga before the release of the laboratory strain (L-Ng19, 5 September), did not cluster with F-Ng21 (the same place, 12 September) and with the nearest Louga H. hebetor but with some from MHM larvae collected in Tondikouarey and Tondibia, located about 12 km from the Louga release point, raising the question of either the representativity of some samples or the real similarity of the sampling conditions. Likewise, all H. hebetor collected in the Dosso site, located in the southern corner of Niger, did not cluster with any others H. hebetor. The level of differentiation from the other sites and the genetic diversity observed point to a low level of gene flow and supply of new alleles to the Dosso subpopulation. The Dosso region is characterized by sedentary livestock breeding and multicrop susbsistence farms no bigger than a few acres, like many areas in Niger, so this sample might be a small isolated one. On the other hand, this finding again underlines the need for more complete and reliable information about each sample to be able to assess them correctly.
To conclude, at both regional and local scales, H. hebetor appeared to be structured in distinct genetic subpopulations, i.e. patches between which gene flow and specimen dispersion may occur by human-mediated and natural short-distance dispersal.
The Bracon brevicornis taxonomic confusion
Phylogenetic reconstructions based on mitochondrial markers and previously conducted cross-species amplifications with the microsatellite loci used in the present study (Garba et al., Reference Garba, Loiseau, Benoit and Gauthier2016), clearly suggest that our putative « B. brevicornis » sample may actually be an H. hebetor one. Indeed, a low level of divergence was identified between the mitotypes obtained for the two species (1–1.3% for COI, 0.3–2.3% for 16S rDNA) leading the B. brevicornis sequences to cluster with H. hebetor sequences. Such a low level of divergence (0.16%) between mtCOI sequences derived from this braconid and low bootstrap values at the taxonomic level making it possible to distinguish between them has already been reported (Rukhsana & Sebastian, Reference Rukhsana and Sebastian2015). In fact, morphologically similar, both using the same wide range of moth hosts, and possibly successfully mating with H. hebetor, B. brevicornis was long considered as a synonym of H. hebetor (Puttarudriah & Basavanna, Reference Puttarudriah and Channa Basavanna1956; Papp, Reference Papp2012). At this time, we can only speculate about the identity of the B. brevicornis sample we studied, but this ongoing taxonomic problem should be addressed not only out of taxonomic interest but also in the framework of conservation BC implemented using local parasitoid communities, as already suggested in previous studies (eg. Sow et al., Reference Sow, Brévault, Delvare, Haran, Benoit, Cœur d'Acier, Galan, Thiaw, Soti and Sembène2018).
Implications for H. hebetor as a biological control agent – conclusions
The reproductive system, based on single-locus sex determination (sl-CSD) of H. hebetor appears to be the driver of genetic diversity and outbreeding of the species in the field. In laboratory strains, the effects of sl-CSD on genetic diversity may be counterbalanced by inbreeding and bottlenecks which, in turn, may lead to a more male-biased sex ratio and declining laboratory rearings. Many studies have demonstrated that sl-CSD substantially influences the sex ratio and the population growth of parasitoid species. Both factors affect the BC potential of parasitoids and as a result, sl-CSD could seriously impede BC if appropriate precautions are not taken (van Wilgenburg et al., Reference Van Wilgenburg, Driessen and Beukeboom2006; Heimpel & de Boer, Reference Heimpel and de Boer2008; Elias et al., Reference Elias, Mazzi and Dorn2009). Despite these considerations, H. hebetor already appears to be able to establish itself when released and is capable of efficient levels of pest control, in particular of MHM in millet crops (Payne et al., Reference Payne, Tapsoba, Baoua, Ba, N'Diaye and Dabiré-Binso2011; Baoua et al., Reference Baoua, Ba, Amadou, Kaboré and Dabire-Binso2018; Kaboré et al., Reference Kaboré, Ba, Ba, Dabire-Binso and Sanon2017). Our study showed that only rare H. hebetor samples collected in the field displayed signatures of a bottleneck, meaning either they quickly re-expanded and grew, which is consistent with their short generation time under Sahel climatic conditions, or they were not seriously affected by the long-lasting diapause of MHM. As a parasitoid of many lepidopteran moths and with the absence of proven diapause, H. hebetor may develop on alternative host species both in the field and in grain stores to survive the period of MHM scarcity. The Louga study region suggests possible continuous connectivity between hosts, and also shows how an area with diverse vegetation may favour H. hebetor genetic diversity and efficiency. Outbreeding and both human-mediated and natural dispersal probably help maintain allelic genetic diversity at the sex locus. Our findings also reveal a complex population structure, probably in patches mainly based on genetically distinct individuals, which, in turn, points to more sexual outcrossing than we might have expected. The possible advantage of keeping a large number of distinct patches – even if they were shown to exchange genes – would be that if many alleles were lost in each patch, at least two different sex alleles would be retained. If so, breeding small numbers of reared individuals rather than larger numbers could be a more efficient way to maintain genetic diversity at the sex-determination locus in H. hebetor. Several times, the question of the natural capacity of H. hebetor to disperse was raised. Our results suggest that this species is weakly dispersed, in any case less than described in the literature but this point needs to be clarified. Neither haplotypes nor alleles typical of laboratory H. hebetor strains were found, which would enable us to adequately monitor the H. hebetor released and their dispersal in the field as well as their dynamics of parasitism on MHM or other hosts. Indeed, even though the releases of H. hebetor help increase parasitism, the observed parasitism cannot only be attributed to augmentative releases. H. hebetor is endemic in most African countries and increases in parasitism result from the cumulative effect of natural and released H. hebetor in the agroecosystems concerned. But, based on microsatellite loci and using a ‘pure’ laboratory strain reared for many generations with no added field specimens, it might be possible to discriminate between the offspring of a released strain and a field strain. Even if some points need to be completed, taken together, our results advance our knowledge of the patterns as well as of the main driver(s) of population genetic variability and structure, demographic events and gene flow of H. hebetor, a widely used BC agent. This study underlines the need to bear the species’ mating system in mind and to be aware of the history of the sample to make sense of the results; it also provides a framework for further population genetics and dynamics studies of BC agents exhibiting CSD such as H. hebetor.
Supplementary material
The supplementary material for this article can be found at https://doi.org/10.1017/S0007485319000142.
Acknowledgements
This study was partly funded by the West Africa Agricultural Productivity Program (WAAP/PPAAO) in Niger and the Institut de Recherche pour le Développement (IRD), France through a collaborative programme between ‘la Direction Générale de la Protection des Végétaux’, Niamey and the IRD. The authors thank the technicians who helped M. Garba to collect the H. hebetor in Niger and E. Rakhshani (University of Zabol, Iran) who graciously sent individuals. The authors are particularly grateful to G. Delvare (CBGP-CIRAD, Montferrier-sur-Lez, France) who helped us check the specific identification of some of the specimens collected and F. Bénédet (UR ‘Forêts et Sociétés’, CIRAD, Montferrier-sur-Lez, France) who helped make the maps. The data used in this work were produced thanks to the technical facilities of the ‘Centre Méditerranéen Environnement Biodiversité’ in Montpellier, France.
Disclosure
The authors declare that they have no conflict of interest.
