Genetic Structure of the Aphis craccivora (Hemiptera: Aphididae) From Thailand Inferred From Mitochondrial COI Gene Sequence.

The cowpea aphid, Aphis craccivora Koch (Hemiptera: Aphididae), is one of the most destructive insect pests of legume plants worldwide. Although outbreaks of this pest occur annually in Thailand causing heavy damage, its genetic structure and demographic history are poorly understood. In order to determine genetic structure and genetic relationship of the geographic populations of this species, we examined sequences of mitochondrial cytochrome c oxidase subunit I (COI) gene of 51 individuals collected from 32 localities throughout Thailand. Within the sequences of these geographic populations, 32 polymorphic sites defined 17 haplotypes, ranging in sequence divergence from 0.2% (1 nucleotide) to 2.7% (16 nucleotides). A relatively high haplotype diversity but low nucleotide diversity was detected in the populations of A. craccivora, a finding that is typical for migratory species. Phylogenetic analysis revealed a weak phylogeographic structuring among the geographic populations and among the haplotypes, indicating their close relationship. Considering the distance between the sampling sites, the occurrence of identical haplotypes over wide areas is noteworthy. Moreover, the low genetic distance (FST ranging from -0.0460 to 0.3263) and high rate of per-generation female migration (Nm ranging from 1.0323 to 20.3333) suggested population exchange and gene flow between the A. craccivora populations in Thailand.

Cowpea aphid, Aphis craccivora Koch, a serious insect pest in legume agriculture, is reported on all continents except the Antarctic (Pettersson et al., 1998, Obopile and Ositile, 2010, Kamphuis et al., 2012). This aphid species is a pest of a wide range of legume species, including yardlong bean (Vigna unguiculata [L.] Walp. ssp. sesquipedalis [L.] H. Ohashi.), cowpea (V. unguiculata [L.] Walp.], chickpea (Cicer arietinum L.), peanuts (Arachis hypogaea L.), garden pea (Pisum sativum L.), and winged bean (Psophocarpus tetragonolobus DC.) (Singh and Jackai, 1985, Annan et al., 1997, Herselman et al., 2004), as well as on members of at least 19 other plant families (Blackman and Eastop, 2007). Alate adults (winged females) are produced whenever food is in short supply and population density is high, as well as during changes in temperature (Blackman and Eastop, 2000). The fecundity, growth, and longevity of A. craccivora vary with physical factors such as weather conditions, moisture, and host plant (Ofuya, 1997).

Significant damages and consequently economic losses to legume crops by A. craccivora are caused by both nymphs and adults (Ofuya, 1997, Blackman and Eastop, 2000) either directly by feeding or indirectly through transmission of at least 14 plant viral diseases, such as cowpea aphid-borne mosaic virus, groundnut rosette virus, peanut mottle virus, peanut stunt virus, and bean common mosaic virus (Sainsbury et al., 2010, Kamphuis et al., 2012). During the seedling stage, the aphids feed on stems and terminal shoots, and later as the plants mature they move to flowers and pods (Ofuya, 1997, Obopile and Ositile, 2010). High infestation with cowpea aphids results in stunted plant growth and delayed flower initiation (Blackman and Eastop, 2000), and if the infestation is not controlled, around a half of crop yield may be lost (Obopile, 2006).

The genetic aspect of the geographical variation in aphid populations can provide important biological information for deploying aphid-resistant cultivars and developing convenient management methods for pest control (Xu et al., 2011). Most aphids exhibit a complicated lifecycle and adapt to diverse host plants (Haack et al., 2000, Ruiz-Montoya et al., 2003), implying that any given geographical population of aphids will be heterogeneous (Lushai and Loxdale, 2002). In addition, different aphid species face different pressure that leads to migration and thus adopt a different migration behavior (Llewellyn et al., 2003, Loxdale and Lushai, 2007, Loxdale et al., 2017). Long-distance migration in aphids depends on wind (Venette and Ragsdale, 2004) because of the poor flight abilities of winged forms (Zhang et al., 2003). Additionally, clonal selection is one of the most important factors in the study on population genetic structure of aphids (Ruiz-Montoya et al., 2003). In aphids, ecological adaption may occur rapidly due to the propagation of asexual parthenogenetic offspring (Lushai and Loxdale, 2003). There are many factors affecting clonal selection of aphids, including host plant, natural enemy pressure, climatic, microclimates, crop density, and resistance to pesticides (Vorburger, 2006, Margaritopoulos et al., 2009, Xin et al., 2014). Moreover, the methods used to control aphid populations, such as pesticide application, can lead to new geno-phenotypes (Xin et al., 2014, Loxdale et al., 2017). Quantification of the genetic diversity of the populations is therefore needed for the management of aphids (Marcelo and Vieira, 2010).

Mitochondrial DNA (mtDNA) has a high mutation rate relative to that of the nuclear genome (Bae et al., 2001, Xu et al., 2012). In addition, it is maternally inherited, does not undergo genetic recombination, and many efficient PCR primers are available (Simon et al., 1999, Hebert et al., 2004, Xu et al., 2012). These properties render the mtDNA valuable for determining genetic structure, gene flow, and evolutionary history of animals (Wilson et al., 1985, Xu et al., 2011). Therefore, such gene sequences have often been used to examine aphid phylogeny (Kim and Lee, 2008), species identification (Foottit et al., 2008, Kim et al., 2010, Shufran and Puterka, 2011), and even genetic structure of aphid populations (Llewellyn et al., 2003, Ruiz-Montoya et al., 2003, Figueroa et al., 2005, Xu et al., 2011, Kharrat et al., 2014).

Among the several genes in mitochondrial genome, cytochrome c oxidase subunit I (COI) gene has been reported to be appropriate for intraspecific variation analysis due to its high degree of variation (Simon et al., 1994, Hu et al., 2008). Therefore, this study used this molecular marker with the aim of characterizing the potential gene flow, genetic diversity, genetic differentiation, and genetic structure of cowpea aphid populations from different geographical regions in Thailand. The results will obtain provide a foundation for optimizing integrated pest management (IPM) programs by adapting control strategies to biological traits shared by several populations of the same genotype.

Materials and Methods

Aphid Samples

Wingless adult A. craccivora were collected in the field from 32 locations throughout its distribution range in Thailand (Supp Table 1 [online only] and Fig. 1). Fifty-one samples were obtained from colonies living on four common host plants, including yardlong bean (Vigna unguiculata subsp. sesquipedalis), peanut (Arachis hypogaea), soybean (Glycine max [L.] Merr.), and winged bean (Psophocarpus tetragonolobus). At least 30 individuals were collected from each collection site. To minimize the chance of sampling individuals from the same colony and hence genotype (Loxdale, 2009, Martens et al., 2009), at each collection site, aphid specimens were taken from a single host plant separated by at least 20 m from the next sample (Xu et al., 2011). Aphids that were collected from the same location were considered one population (Xin et al., 2014). All A. craccivora samples were preserved in 100% ethanol and kept at −20 °C until DNA extraction.

Table 1.

The 28 polymorphic nucleotide sites among 17 mitochondrial cytochrome c oxidase I (COI) haplotypes of Aphis craccivora in Thailand

	58	99	138	147	276	313	348	361	405	420	425	465	493	503	506	507	510	511	512	513	516	524	541	553	556	558	559	567	582	586	596	597	599
HAP1	G	A	C	A	A	G	T	T	G	A	A	A	G	T	C	T	A	T	T	A	A	G	T	G	G	G	G	T	A	T	T	T	T
HAP2	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	T	.	T	.	.	.	.	G	.	.	C	.
HAP3	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	T	.	.	.	.	.	.	.	.	.
HAP4	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	A	T	A	.	.	.	G	.	.	.	A
HAP5	.	.	.	.	.	.	.	.	.	.	.	.	T	A	G	A	T	A	A	T	.	.	.	T	C	A	T	.	.	.	C	C	A
HAP6	.	.	.	.	.	.	.	C	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	A
HAP7	A	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.
HAP8	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	G	.	.	T	.	.	.	.	.	A	.	C	A
HAP9	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	T	.	.	.	.	.	.	.	.	A
HAP10	.	G	.	.	.	A	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	T	.	.	.	.	.	.	.	.	.
HAP11	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	T	C	A	T	.	.	.	C	C	A
HAP12	.	.	.	.	.	.	.	C	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.
HAP13	.	.	T	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.
HAP14	.	.	.	.	.	.	.	C	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	T	.	.	.	A	.	.	.	.	A
HAP15	.	.	.	.	.	.	G	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.	.
HAP16	.	.	.	.	.	.	.	.	.	.	G	.	.	.	G	A	.	A	.	T	.	.	.	T	C	A	T	.	.	.	C	C	A
HAP17	.	.	.	T	T	.	.	.	T	T	.	T	.	.	.	.	.	.	.	T	.	.	.	T	C	A	T	.	.	.	C	C	A

Only positions that are different from haplotype HAP1 are indicated.

Fig. 1.

Location of the sampling sites of Aphis craccivora in Thailand. The numbers of the samples correspond to those in Supp Table 1 [online only].

DNA Extraction, Amplification, and Sequencing

Total genomic DNA was extracted from wingless adult aphid individuals using a DNeasy Blood & Tissue kit (Qiagen, Germantown, MD) according to the manufacturer’s instructions. A part of the protein-coding mitochondrial COI gene was amplified by PCR using primers LCO1490 (5ʹ-GGTCAACAAATCATAAAGATATTGG-3ʹ) and HCO2198 (5ʹ-TAAACTTCAGGGTGACCAAAAAATCA-3ʹ) (Folmer et al., 1994). PCRs were performed in a final volume of 40 µl containing 1× PCR master mix (Fermentas Thermo Fisher Scientific, Waltham, MA), 20 µmol/l of each primer, and at least 20 ng of genomic DNA template. The thermocycling profiles consisted of an initial denaturation step at 94 °C for 5 min, followed by 35 cycles of denaturation at 94 °C for 1 min, annealing at 48 °C for 1 min; and extension at 72 °C for 150 s with a final extension at 72 °C for 5 min. The PCR products were visualized by electrophoresis on a 1.5% agarose gel, and those with a single band were purified by using QIAquick Gel Extraction kit (Qiagen). The purified PCR products were sequenced directly at AITbiotech Pte Ltd (Singapore).

Sequence Alignment and Diversity Indices

The raw DNA sequences were each checked visually and verified for protein coding frame-shifts to avoid pseudogenes (Zhang and Hewitt, 1996). The sequences obtained in the present study were deposited in the National Center for Biotechnology Information (NCBI) GenBank, and the accession numbers are listed in Supp Table 1 [online only]. MEGA6 v6.06 (Tamura et al., 2013) was used to generate multiple alignments and to calculate the number of substitutions and nucleotide compositions for COI. The genetic variation was then investigated as the average number of nucleotide differences (k), number of polymorphic (segregating) sites (S), number of haplotypes (No.), haplotype diversity (hd) (Nei, 1987), and nucleotide diversity (Pi) (Nei and Li, 1979) by using DNAsp v5.0 program (Librado and Rozas, 2009).

Historical Demography

To examined historical demography of A. craccivora in Thailand, we calculated Tajima’s D (Tajima, 1989) and Fu’s Fs (Fu, 1997) neutrality test as an assessment of possible population expansion as a deviation from the neutrality. Both statistics were performed using ARLEQUIN ver 3.5.2.2 (Excoffier and Lischer, 2010) and significance was assessed with 1,000 bootstrap replicates. A significantly positive value of Tajima’s D indicates population subdivision or population contraction, while a significantly negative value indicates recent population size expansion (e.g., after a bottleneck). Regarding the Fu’s Fs value, a large negative value indicates the excess of rare alleles in the population, which implies a recent increase in population size. Additionally, the Ramos-Onsins and Rozas’s R2 (Ramos-Onsins and Rozas, 2002) was performed using DNAsp v5.0 program (Librado and Rozas, 2009) and significance was assessed with 1,000 coalescent simulations. All statistical analyses were performed for each geographical area separately and for each haplotype.

Phylogenetic Analysis

For phylogenetic analysis, we used maximum likelihood (ML) and Bayesian inference (BI) methods to reconstruct the phylogenetic relationships among COI haplotypes and among geographical samples. Prior to ML and BI analyses, the program Kakusan4 (Tanabe, 2007) was used to determine the best-fit models of nucleotide substitution as judged by the Akaike information criterion (AIC) (Akaike, 1974) for ML and the Bayesian information criterion (BIC) (Schwarz, 1978) for BI. The ML analysis was performed in Treefinder (Jobb et al., 2004), using likelihood-ratchet method with 1,000 bootstrap replicates to examine branch confidence values. The topology that received bootstrap support (BS) of at least 70% was considered sufficiently resolved (Huelsenbeck and Hillis, 1993). The BI analysis was performed using MrBayes v3.1.2 (Huelsenbeck and Ronquist, 2001). The Bayesian analysis using the selected model consisted of four simultaneous chains run for 3,000,000 generations starting from a random tree, with tree sampling every 100 generations.

Genetic Distance and Migration Estimate

The pairwise genetic distance (FST) and per-generation female migration rate (Nm) of Thai A. craccivora populations were estimated using ARLEQUIN ver 3.5.2.2 (Excoffier and Lischer, 2010). Population pairwise genetic distance was estimated with 1,000 bootstraps following Excoffier et al. (1992). Moreover, the distances between DNA sequences were calculated by Kimura’s 2-parameter method (Kimura, 1980). The pairwise FST values were used to examine per-generation female migration rate as follow: Nm = (1−FST)/4FST.

Mantel Test for Isolation by Distance

To examine the relationship between genetic and geographical distances, a Mantel’s test (Mantel, 1967) with 1,000 permutation was implemented using the program IBD version 1.5.2 (Bohonak, 2002). Population pairwise FST was used for genetic distance, and geographical distance between populations were calculated according to the coordinated (latitude and longitude) of the location of each population (Hufbauer et al., 2004, Li et al., 2015). The natural logarithm of distance was used as the measure of geographical distance to minimize error (Slatkin, 1993, Sokal and Rohlf, 1995).

Results

COI Gene Sequence Analysis

After removing the primers used in the PCRs, the effective COI sequences consisted of 602 nucleotides. The analysis of the sequences revealed a total of 17 haplotypes (designated HAP1–HAP17) in samples collected at the 32 localities throughout Thailand. Of the 51 individual A. craccivora specimens, 26 individuals were included in HAP1 and 8 individuals possessed HAP12, together constituting 66.67% of aphid specimens (Supp Table 1 [online only]). These two haplotypes were found in all locations studied. All 17 haplotypes showed 33 polymorphic sites, of which 7 were A/G transitions, 4 T/C transitions, 14 A/T transversions, 6 G/T transversions, 1 G/C transversion, and 1 transition and transversion (Table 1). The percentage of parsimony-informative characters relative to the total number of characters for COI gene was 2.32%, and the frequency of transition/transversion mutations was 36.36/63.64. Nucleotide composition of COI was 35.3% A, 40.6% T, 10.2% G, and 13.8% C. Therefore, the A + T nucleotide bias reported in the animal mitochondrial genome (Boore, 1999) was also found in the A. craccivora COI sequences, as also detected from previous studies on aphids (Simon et al., 1994, Von Dohlen and Moran, 2000, Von Dohlen et al., 2002, Lee et al., 2011, Lassaad et al., 2013).

Population Genetic Structure and Haplotype Divergence

The summary of molecular diversity indices of A. craccivora from Thailand are given in Table 2. Of the 17 haplotypes, 13 were unique to one individual and three were shared by individuals from different populations. The haplotype diversity (hd), nucleotide diversity (Pi), and average number of nucleotide differences (k) of all samples were 0.721 ± 0.064, 0.00511 ± 0.00129, and 2.676, respectively. Based on mitochondrial COI gene sequences, all populations revealed medium to high haplotype diversity (ranging from 0.154 to 1.000) but relatively low nucleotide diversity (ranging from 0.00102 to 0.01063) (Table 2). The samples from the central regions possessed the highest haplotype and nucleotide diversity, while samples from the northern part had the lowest haplotype (0.154 ± 0.126) and nucleotide diversity (0.00102 ± 0.00084).

Table 2.

Summary of molecular diversity indices and population expansion test statistics of mitochondrial cytochrome c oxidase I (COI) genes

Localities	N	No.	S	k	hd (±SD)	Pi (±SD)	D	Fs	R2
North	13	2	4	0.615	0.154 (0.126)	0.00102 (0.00084)	−1.775*	1.474	0.267
Central	10	10	23	6.400	1.000 (0.045)	0.01063 (0.00268)	−1.167	−4.947*	0.125
Northeast	9	3	5	1.667	0.639 (0.126)	0.00277 (0.00086)	−0.398	1.378	0.213
West	9	4	14	3.278	0.694 (0.147)	0.00544 (0.00325)	−1.746*	1.505	0.268
South	10	2	13	2.60	0.200 (0.154)	0.00432 (0.00333)	−1.976**	4.861	0.300
All samples	51	17	33	2.676	0.721 (0.064)	0.00511 (0.00129)	−1.987**	−5.837*	0.046*

Number of individuals (N), number of haplotypes (No.), number of polymorphic (segregation) sites (S), average number of nucleotide differences (k), haplotype diversity (hd) and nucleotide diversity (P) with standard deviation (SD), Tajima’s D, Fu’s F and Ramos-Onsins and Rozas’ R2. * and ** indicate significant difference at P<0.05 and P<0.01, respectively.

The pairwise distance comparisons among A. craccivora haplotypes are shown in Table 3. Sequence divergence among the 17 haplotypes by pairwise comparison ranged from 0.2% to 2.7% (ranging from 1 to 16 bp). The highest sequence divergence was found when haplotype HAP5 was compared with HAP10, HAP12, HAP13, and HAP15 (Table 3). When compared with other aphid species, very high sequence divergence of 7.0–9.8% (ranging from 40 to 55 bp) was detected among the 17 COI haplotypes of A. craccivora and A. glycines (AB560725) (Table 3). These results confirmed that intraspecific sequence divergence was lower than interspecific sequence divergence.

Table 3.

Pairwise comparison of nucleotide sequences of the partial mitochondrial cytochrome c oxidase I (COI) gene of Aphis craccivora populations in Thailand

	1	2	3	4	5	6	7	8	9	10	11	12	13	14	15	16	17	18
1. HAP1	–	0.007	0.002	0.008	0.025	0.003	0.002	0.008	0.003	0.005	0.012	0.002	0.002	0.007	0.002	0.020	0.022	0.071
2. HAP2	4	–	0.005	0.008	0.025	0.010	0.008	0.008	0.007	0.008	0.012	0.008	0.008	0.010	0.008	0.020	0.022	0.079
3. HAP3	1	3	–	0.007	0.024	0.005	0.003	0.007	0.002	0.003	0.010	0.003	0.003	0.005	0.003	0.019	0.020	0.073
4. HAP4	5	5	4	–	0.025	0.008	0.010	0.010	0.005	0.010	0.012	0.010	0.010	0.008	0.010	0.020	0.022	0.081
5. HAP5	15	15	14	15	–	0.025	0.027	0.024	0.022	0.027	0.013	0.027	0.027	0.025	0.027	0.008	0.020	0.098
6. HAP6	2	6	3	5	15	–	0.005	0.008	0.003	0.008	0.012	0.002	0.005	0.003	0.005	0.020	0.022	0.075
7. HAP7	1	5	2	6	16	3	–	0.010	0.005	0.007	0.013	0.003	0.003	0.008	0.003	0.022	0.024	0.073
8. HAP8	5	5	4	6	14	5	6	–	0.005	0.010	0.010	0.010	0.010	0.008	0.010	0.019	0.020	0.081
9. HAP9	2	4	1	3	13	2	3	3	–	0.005	0.008	0.005	0.005	0.003	0.005	0.017	0.019	0.075
10. HAP10	3	5	2	6	16	5	4	6	3	–	0.013	0.007	0.007	0.008	0.007	0.022	0.024	0.077
11. HAP11	7	7	6	7	8	7	8	6	5	8	–	0.013	0.013	0.012	0.013	0.008	0.010	0.083
12. HAP12	1	5	2	6	16	1	2	6	3	4	8	–	0.003	0.005	0.003	0.022	0.024	0.073
13. HAP13	1	5	2	6	16	3	2	6	3	4	8	2	–	0.008	0.003	0.022	0.024	0.070
14. HAP14	4	6	3	5	15	2	5	5	2	5	7	3	5	–	0.008	0.020	0.022	0.079
15. HAP15	1	5	2	6	16	3	2	6	3	4	8	2	2	5	–	0.022	0.024	0.073
16. HAP16	12	12	11	12	5	12	13	11	10	13	5	13	13	12	13	–	0.015	0.092
17. HAP17	13	13	12	13	12	13	14	12	11	14	6	14	14	13	14	9	–	0.092
18. A. glycines	41	45	42	46	55	43	42	46	43	44	47	42	40	45	42	52	52	–

Numbers above the diagonal are mean distance values; numbers below the diagonal are absolute distance values. GenBank accession number of Aphis glycines is AB506725 (Komazaki et al., 2010).

Demographic History

Demographic history of A. craccivora in Thailand is demonstrated in Table 2. When considering the A. craccivora samples in Thailand as a single population, the value of Ramos-Onsins and Rozas’ R2 was small and positive (R2 = 0.046; P = 0.011), indicating A. craccivora population growth. Additionally, Tajima’s D was significantly negative (D = −1.987; P = 0.002), revealing a recent population size expansion (e.g., after a bottleneck or a selective sweep). In addition, the large significantly negative value of Fu’s Fs (Fs = −5.837; P = 0.019) indicated the excess of rare alleles in the population of A. craccivora of Thailand. This implies a recent increase in population size or positive selection (Table 2). When separate the samples into five geographical populations, the significantly negative Fu’s Fs value (Fs = −4.947; P = 0.005) suggests a recent population expansion in the central region. The values of Tajima’s D were negative and those of Ramos-Onsins and Rozas’ R2 positive, but none of them were significant. In summary, most A. craccivora populations in Thailand revealed signs of expansion, although not always at a significant level.

The best-fit model for ML analysis under AIC was TVM + Gamma, while the best-fit model of evolution for BI analysis under BIC was HKY85 + Homogeneous. Gene genealogies estimated by ML and BI were highly congruent, and therefore, only the tree topologies from ML analyses are presented herein. Phylogenetic relationships among COI haplotypes and among geographical regions of A. craccivora are shown in Fig. 2A and B, respectively. The phylogenetic tree from the mitochondrial COI gene strongly supported the monophyly of A. craccivora (BS = 100%, Bayesian posterior probability [PP] = 1.0). The 17 haplotypes were resolved into five clusters, most of which were weakly supported (BS ranging from 52% to 67%; PP ranging from 0.63 to 0.67). Such a low support might be due to small nucleotide differences among the haplotypes (Tables 1 and 3).

Fig. 2.

Phylogenetic relationships among cytochrome c oxidase I (COI) haplotypes (A) and among geographical regions (B) of Aphis craccivora in Thailand based on maximum likelihood (ML). Aphis glycines (AB506725) was utilized as outgroup. Node supports inferred from Bayesian posterior probability and bootstrap value for ML.

The phylogenetic analysis of the five geographical regions resolved four clades (Fig. 2B), most of which were weakly supported (BS < 60%). Clade I includes the aphid samples from the northeast, central, and west, Clade II includes most of the samples from the north, west, and south regions, Clade III consists of aphid samples from all the regions, and Clade IV comprises the samples from the northeast regions. However, there was no evidence for strong geographical clustering in the phylogenetic tree (Fig. 2B).

The pairwise genetic distance (FST) and per-generation migration rate (Nm) of A. craccivora populations in Thailand are shown in Table 4. The pairwise FST ranged from −0.0460 to 0.3263; only six pairs of populations showed non–statistically significant genetic differentiation (P > 0.05). The pairwise per-generation female migration rate (Nm) between all populations was greater than one (Table 4), indicating that more than one female A. craccivora per generation was estimated to migrate between all population pairs and showing that there is at least some ongoing gene flow between most A. craccivora populations in Thailand.

Table 4.

Genetic distance (FST) and gene flow (Nm) among populations of Aphis craccivora

	North	Central	Northeast	West	South
North	–	0.2254**	0.3263*	0.0240	−0.0460
Central	1.1788	–	0.1801*	0.0550	0.0982
Northeast	1.0323	2.2762	–	0.0686	0.1723**
West	20.3333	8.5091	6.7886	–	−0.0526
South	11.3695	4.5916	2.4019	10.0057	–

Inf, infinite.

Values above the diagonal are the estimates of genetic distance (FST) and values below the diagonal are the estimates of per-generation female migration rate (Nm).

and **indicate significant difference at P < 0.05 and P < 0.01, respectively.

Isolation by Distance Among Populations

In this study, geographical distance between populations was calculated according to the latitude and longitude of each population. The geographical distance between pair of populations ranged from 5.86 to 1,279 km. For the relationship between pairwise FST and geographic distance, a significant positive pattern of isolation–by–distance was found among the Thai’s aphid populations (Fig. 3; Z= 91.448, r = 0.2122, P = 0.017).

Fig. 3.

Pairwise genetic distance (FST) plotted against the log-transformed geographic distances (in kilometer) between Thai’s aphid populations. The solid line shows the best-fit linear regression based on all points.

Discussion

In the present study, the maximum pairwise divergence of COI gene sequences within A. craccivora was 2.7%. In previous studies that utilized the partial COI gene sequences in aphids, this divergence was 1.7% for Sitobion avenae (Fabricius) (Xu et al., 2011) and 0.4% for Acyrthosiphon pisum (Harris) (Boulding, 1998). Even most aphids are plant specialists, but by using DNA markers, the significant genetic structure differences can be expected including inter- and intraspecific variation (Loxdale, 2009), cryptic (Sunnucks et al., 1997) hybridization and interclonal, intermorph genetic differences (Lushai et al., 1997, Loxdale and Harvey, 2016). Compared to other reported values, the magnitude of sequence divergence in A. craccivora which is therefore higher than that in other aphid species, possibly revealed cryptic hybridization in this species.

The genetic diversity of aphids is influenced by several factors, such as host plant, topography, and climate (Cai et al., 2004). Our study revealed a relatively high haplotype diversity (hd) and low nucleotide diversity (Pi) in the populations of A. craccivora from Thailand (Table 2), which might be exhibitive of genetic bottleneck events (Wei et al., 2013). In particular, the haplotype diversity of the central population was higher than the diversity of other populations. This may be due to the greater gene flow between this geographic location and the surrounding legume plantations.

Phylogenetic analysis of A. craccivora identified no haplotype group or a notable polymorphic haplotype in Thai populations. Low node support in both topologies based on geographical populations and haplotypes (Fig. 2) indicated close relationship among the populations. In particular, haplotypes HAP1 and HAP12 were presented at all localities surveyed (Supp Table 1 [online only]). Given the geographical distance between the sampling sites, the occurrence of identical haplotypes over wide areas is noteworthy. However, there was an obscure relationship between genetic and geographic distance. A Mantel test showed a significant positive correlation between genetic and geographical distance (r = 0.2122, P = 0.017; Fig. 3). Thus, we cannot ignore that the geographic isolation is one of importance factors contributing to genetic differentiation of this aphid species. Similar to our result, Cao et al. (2012) found significant isolation–by–distance (r = 0.6392, P < 0.05) in A. spiraecola populations from China. A significant positive relationships between genetic and geographical distance was also detected among cotton aphid (A. gossypii) based on both microsatellite DNA markers (r = 0.433, P = 0.00) and combined mitochondrial COI and Cytb genes sequences (r = 0.150, P < 0.013) (Wang et al., 2017). Frequently, isolation–by–distance effects are pronounced in moderately mobile species but weak in low and high mobility species (Wang et al., 2017). Moderate flight and migration ability of aphids might allow short − distance movement. Therefore, this might be the reason that there is a significant isolation–by–distance among A. craccivora population.

The population analysis showed that most of A. craccivora populations were inter-population genetically differences to each other, suggesting a moderate to high rate of gene flow, certainly enough to counteract the dual effects of selection and drift locally (Table 4). Llewellyn et al. (2003) reported that the interaction of genetic drift, gene flow, migration, natural selection, and mutation affect the genetic structure of populations. No genetic differentiation among populations is for FST = 0, while the populations are completely differentiated if FST = 1. Wright (1978) divided the threshold of FST into four levels as follow: 1) FST > 0.25 indicates significant divergence among populations; 2) 0.15 < FST ≤ 0.25 high genetic differentiation; 3) 0.05 < FST ≤ 0.15 moderate genetic differentiation; and 4) FST ≤ 0.05 low genetic differentiation. Our results showed low to high genetic differentiation of the five geographical populations; significant divergence was detected only between the north and northeast populations (Table 4). The sampling sites in the present study covered a wide range of major legume plantations (32 localities) in Thailand. Although the straight line distance between Chiang Rai (north) and Krabi (south) is ∼1,600 km, most of the population pairs showed no significant difference in genetic distance (FST), arguing that most of A. craccivora populations formed a single genetic group.

Our results also indicated that A. craccivora genotypes on the less suitable host plants in surrounding area acting as a reservoirs for aphid population did not significantly contribute to population build-up in the legume fields during the growing season, but the migrant aphid populations are greater considerable for the establishment of the population (Klueken et al., 2012, Xin et al., 2014). Thus, the long distance migration of aphids might be an important factor in population genetic structure. Aphids can be suspended in the air and migrate long distances with the wind (Delmotte et al., 2002), but are always actively flying even so, as the late, great L.R. (Roy) Taylor showed many moons ago (Loxdale et al., 1993).

If the migration distance of the aphids is large enough, the flight phenomenon can overcome the forces of genetic drift and mutation (Llewellyn et al., 2003, Xin et al., 2014). It is known that when the high population density and host plant cannot supply adequate nutrition, alate (winged) aphids are produced, facilitating long distance migration (Blackman and Eastop, 2000). Long distance migration in aphids was reported in several publications (Loxdale and Lushai, 2007). The molecular technique is evidently used for long distance migration in aphids eg. Sitobion avenae (Fabricius) (Llewellyn et al., 2003). However, the migration of A. craccivora in relation to molecular marker is still under-research. In addition, the per-generation female migration (Nm) was greater than one in all pairwise comparisons (ranging from 1.0323 to 20.3333); high Nm of 11.3695 was detected between the farthest distance (north and south) (Table 4). According to Wright (1951), the Nm values >1.0 indicate a low genetic differentiation among populations, whereas those <1.0 the opposite. Thus, our results suggest that dispersal over long distance and with no geographical barriers are major factors in the demography of A. craccivora.

Interestingly, high genetic differentiation was detected between the north and northeast populations, which might result from a limited dispersal of A. craccivora between the two populations. The gene flow between aphid populations in the north and northeast regions might be impacted by the barrier created by the Luang Prabang Range Mountains, Doi Khun Tan Mountains, and Phetchabun Mountains, all of which separate northern, central, and northeastern regions of Thailand. This phenomena of genetic isolation by geographical barriers were also reported in Sitobion avenae (Fabricius) (Xin et al., 2014), Diuraphis noxia (Kurdjumov) (Zhang et al., 2012), Macrosiphoniella tanacetaria (Kaltenbach) and Metopeurum fuscoviride (Stroyan) (Loxdale et al., 2011).

The A. craccivora has been considered as one of the most important migratory pest insects that is widely distributed in Thailand throughout the regions cultivated with different legume species, especially V. unguiculata (L.) Walp. ssp. sesquipedalis (L.) H. Ohashi. and V. radiata (L.) R. Wilczek. Herein, we examined the population genetic diversity, population genetic differentiation, population genetic structure, phylogenetic relationship, historical demography, and genetic distance and estimated migration of five different A. craccivora geographical populations in Thailand based on the sequences of mitochondrial COI. Our results, we hope and trust, provide a basis to improve IPM programs through developing control strategies that target biological traits shared by various populations.

Supplementary Data

Supplementary data are available at Journal of Insect Science online.

Compliance With Ethical Standards

Conflict of interest: The authors have no conflict of interest.

Ethics statement: The cowpea aphid (Aphis craccivora) is one of the most serious insect pests of many legume plant species. It is not endangered species or protected by law. Therefore, no permits are required to study this insect.

Click here for additional data file.