TY - JOUR
AU - Brown, Judith, K
AB - Abstract The association between Bemisia tabaci mitotypes and cotton leaf curl outbreaks in Pakistan was investigated using the mitochondria cytochrome oxidase I gene (COI) as a molecular marker. The 3′-651 base fragment has been used to resolve B. tabaci phylogenies. However, the 5′-618 base fragment was nearly unexplored. Phylogenetic analysis for 829 whiteflies from 11 districts in two provinces of Pakistan, indicated all haplotypes grouped on the Asia II major clade, with Asia II-1 mitotype predominating, at 84%, compared to Asia II-5 and II-7, at ~16%, combined. The 3′- and 5′-fragment tree topologies were similar, while the concatenated topology was unique in some respects. Comparisons of segregating sites within the 3′- and 5′-loci, at third codon positions, 71 and 47, and of transitions to transversions (Ti/Tv) ratio of 2.93 and 5.9, respectively, showed the 3′-locus was most informative, while nucleotide diversity (π) was highest for the 5′-end, indicating both fragments contributed to concatenated tree structure. The extent of haplotype diversity, measured by Tajima’s D, R2, and Fu’s F analyses, revealed significant demographic expansion for Asia II-1 and II-7 mitotypes. The bottleneck that preceded the expansions was evident in the temporal changes in mtCOI polymorphisms beginning in ~1990s, a timeframe known to have coincided with the adoption of a high-yield whitefly-susceptible cultivar in 1988, followed by pesticide overuse. These two cooperating phenomena appear to have exerted selection on the cotton leaf curl disease (CLCuD)-whitefly complex, resulting in the emergence of a resistance-breaking begomovirus as the polyphagous Asia II-1 mitotype underwent a genetic expansion that led to ‘a perfect storm’. Whitefly-transmitted geminiviruses (Begomovirus: Geminiviridae) have emerged as important plant pathogens of cotton in Pakistan, beginning in ~1990, with the outbreak of Cotton leaf curl Multan virus (CLCuMV), a begomovirus causing widespread disease, resulting in 50–80% loss of the crop (Mansoor et al. 2003a). The disease subsided with the development and release of a virus-resistant cotton variety. In 2000–2001, a recombinant, Cotton leaf curl Kokhran virus-Burewala (CLCuKoV-Bu) emerged, causing a second outbreak, apparently driven by the widespread cultivation of the CLCuMV-resistant cotton variety that precipitated the emergence of the resistant-breaking ‘Burewala’ strain (Mahmood et al. 2003, Mansoor et al. 2003a). These two major CLCuD pandemics that have occurred in cotton in Pakistan have been attributed to several factors, including the unprecedented high population sizes of the whitefly vector, and the unexpected rapid diversification of CLCuD-associated begomoviral species and strains (Briddon and Markham 2000). However, a better understanding of the spatio-temporal dynamics of the whitefly vectors that have driven the spread of the CLCuD-begomoviral pathogens is expected to aid in the development of effective management practices. The whitefly Bemisia tabaci (Gennadius) is widely recognized as an important insect pest and virus vector in cotton-vegetable cropping systems (Byrne and Bellows 1991, Brown and Czosnek 2002, Jones 2003) in tropical and sub-tropical world regions (Russell 1957, Mound and Halsey 1978, Rosell et al. 1997). The identification of whiteflies to species is based on a set of defined morphological characters on the fourth instar stadium, or ‘pupa’ (Takahashi 1936, Russell 1948). However, it has been shown that fourth instar nymphs exhibit morphological plasticity, manifest as the development of different numbers and placements of setae in response to different plant leaf surfaces, some of which are diagnostic characters (Mound 1963, Lopez-Avila 1986, Mohanty and Basu 1986), thereby confounding whitefly species identification, and resulting in misclassification. To remedy this problem, three previously recognized genera and 20 species were synonymized, creating the present-day B. tabaci epithet (Russell 1957, Martin and Mound 2007). In contrast to the morphological stasis of this group, biological studies of a set of globally representative B. tabaci has revealed genetic variability with a basis in phylogeographical endemism and the distinct ecological niches within (see references in Brown 2010), as well as differences in host range, insecticide resistance, and plant virus transmission efficiency (Bird 1957; Costa and Russell 1975; Costa and Brown 1991; Costa et al. 1993; Gawel and Barlett 1993; Perring et al. 1993; Brown et al. 1995a; Moya et al. 2001; Legg et al. 2002; Berry et al. 2004; De Barro et al. 2005; Delatte et al. 2005; Sseruwagi et al. 2005, 2006; Brown 2007; Chowda-Reddy et al. 2012; Legg et al. 2014); leading to the classification of the group into ‘biological types’, or ‘biotypes’ (Costa and Brown 1991; Burban et al. 1992; Costa et al. 1993; Bedford et al. 1994; Brown et al. 1995a; Viscarret et al. 2003; Delatte et al. 2005; Sseruwagi et al. 2005, 2006). In addition, the term ‘host race’ has been used to refer to B. tabaci variants to differentiate those with limited host range (Bird 1957, Costa and Russell 1975, Fishpool and Burban 1994, Legg et al. 2002), from others are polyphagous (Bird 1957, Fishpool and Burban 1994). However, in the last two decades, the use of the 3′-fragment of the mitochondria cytochrome oxidase I gene (COI) (~850 bp) as a molecular marker for differentiating whitefly populations worldwide, has allowed the classification of the B. tabaci complex into mitochondrial types (mitotype groups) (Brown et al. 1995b, Frohlich et al. 1999) that group phylogeographically (Brown et al. 1995a,b; De Barro et al. 2000; Legg et al. 2002; Viscarret et al. 2003; Berry et al. 2004; De Barro et al. 2005; Sseruwagi et al. 2005, 2006; Boykin et al. 2007; Brown 2010; Chu et al. 2010; Dinsdale et al. 2010; Gill and Brown 2010; Hu et al. 2011; Alemandri et al. 2012; Chowda-Reddy et al. 2012; Esterhuizen et al. 2013; Firdaus et al. 2013; Lee et al. 2013; Ashfaq et al. 2014). For convention, throughout this article, mitochondrial type (mitotype) is used to refer to a group of closely related haplotypes that cluster in the same lineage within a major phylogeographic clade, and designations of mitotypes at the subgroup level are considered putative, albeit, to enable referencing some recent literature. In contrast, the 5′-COI sequence (Folmer et al. 1994) has almost never been used as a marker for differentiation of the B. tabaci complex, despite its recent popularity for ‘barcoding’ arthropods (Brower 1994, Hebert et al. 2004, Smith et al. 2005, Ward et al. 2005, Hajibabaei et al. 2006, Pons et al. 2006, Ivanova et al. 2007). In one exception, the ‘Folmer fragment’, i.e., the 5′-end of the COI gene was used to characterize B. tabaci variants in Pakistan cotton crops (Ashfaq et al. 2014). Despite the relatively small sample size, the tree topology was similar to those from previous studies for B. tabaci sampled from the same or nearby locales in Pakistan (Simón et al. 2003; Ahmed et al. 2010, 2011), including collections analyzed most recently (Masood et al. 2017, Islam et al. 2018). However, whether differences existed between the phylogenetic signal and substitution patterns for the 5′- and 3′-COI fragments that could influence tree structure and confound the phylogenetic conclusions (Erpenbeck et al. 2006), has not been addressed. Important contributing factors that remain unstudied in the context of the successive leaf curl disease outbreaks in Pakistan, are the prospective influences of population structure and the potential for genetic expansion on the temporal and spatial distributions of extant B. tabaci mitotypes. Based on records that are available from the time of the initial outbreak and spread of CLCuMuV, to the emergence and spread one decade later of CLCuKoV-Bur, and persisting to the present (Mansoor et al. 2003b, Zaidi et al. 2016, Hassan et al. 2017), it is possible to investigate the dynamics of whitefly vector-virus populations. These results could aid greatly in explaining the phylodynamics of the begomovirus emergence in relation to fluctuations in the whitefly mitotypes in the affected areas, in the context of shifts in cotton germplasm genotypes selected to enhance cotton production while also managing pest and disease constraints. Here, the 3′-COI marker, routinely used for haplotyping B. tabaci, and the rarely used ‘Folmer’ fragment located at the 5′-end of the COI gene, were polymerase chain reaction (PCR)-amplified from field-collected B. tabaci, sequenced, and used to investigate the distribution and genetic variability of B. tabaci mitotypes associated with cotton and non-cotton host plants in ten districts where the predominant crop is cotton, or those associated with vegetables-urban landscape near Lahore. Whitefly samples were collected during 2012–2014, a timeframe coinciding with the spread of CLCuKoV-Bu, approximately one decade after this recombinant begomovirus was first discovered near Burewala in the Vehari district of the Punjab Province, Pakistan (Mahmood et al. 2003, Rajagopalan et al. 2012). In this study, the genetic structure based on diversity of the haplotypes was also evaluated temporally, considering characteristics of the B. tabaci populations studied here and from previous investigations to aid in providing insights into the dynamics of mitotypes populations in relation to CLCuD outbreaks in Pakistan. Also, in this study the extent to which the 5′- and 3′-end sequences of the mtCOI gene, respectively, were informative based on polymorphisms, divergence, and tree topology was evaluated, and results were compared for the separate gene fragments and a concatenated sequence consisting of both fragments. Materials and Methods Whitefly Samples and DNA Isolation Whitefly B. tabaci adults were collected from cotton and non-cotton host plants in ten major cotton-growing districts in Pakistan: Sahiwal, Lodhran, Rahim Yar Khan, Bahawalpur, Bahawalnagar, Multan, Vehari, Khanewal, Pakpattan, and Shaheed Benazir Abad, and one non-cotton-growing district: Lahore. Multiple samples were collected from symptomatic cotton, and other cultivated (‘non-cotton’) or uncultivated (wild) host plants showing leaf curling, vein thickening, and/or enations, indicative of a virus-like disease. Each sample consisted of at least six whiteflies per plant. At least three samples were collected per field, and at least three fields were sampled per district. Whiteflies were placed live into a 1.5-ml microfuge tube containing 95% ethanol, and samples were stored at −20°C. Eight hundred twenty-nine adult whiteflies, two–three per field collection, were processed for the study. Total DNA was purified from individual whiteflies using the method of Zhang et al. (1998). Each whitefly was placed onto a piece of filter paper for several seconds to remove excess ethanol, and ground in 600 µl CTAB buffer (100 mM Tris–HCl pH 8.0, 20 mM EDTA pH 8.0, 1.4 M NaCl containing 0.2% 2-mercaptoethanol, 2% hexadecyltrimethyl-ammonium bromide), and incubated at 65°C for 15 min. One vol chloroform was added, mixed, and the emulsion was broken by centrifugation at 12,000 RPM at 4°C for 3 min. The supernatant was transferred to a sterile microfuge tube, and 1 vol isopropanol and 40 μg glycogen were added, with incubation at 4°C for 10 min. The DNA pellet was collected by centrifugation at 12,000 RPM at 4°C for 10 min, washed with 70% ethanol, air-dried, and dissolved in 10 mM Tris buffer, pH 8.0. Amplification, Cloning, and DNA Sequencing The 5′-end spans the nucleotide (nt) coordinates, ~42–659 (618 bases) (Folmer et al. 1994), compared to the 3′-end fragment, which is located at coordinates 699–1450 (651 bases) (Simon et al. 1994, Frohlich et al. 1999). The latter fragments of the whitefly COI gene were amplified by PCR. The 5′-end amplicon of ~700 bp in size, was PCR-amplified using primers: F-LCO1490 5′-GGTCAACAAATCATAAAGATATTGG-3′ and R-HCO2198 5′-TAAAGTTCAGGGTGACCAAAAAATCA-3′ (Folmer et al. 1994), and the 3′-COI fragment, at ~850 bp in size, was amplified using the primers: F-C1-J-2195-5′-TTGATTTTTTGGTCATCCAGAAGT-3′ and R-L2-N-3014 5′-TCCAATGCACTAATCTGCCATATTA-3. This primer pair was designed for the COI amplification of selected insect groups (Simon et al. 1994), a region that was first demonstrated to be an informative molecular marker for the whitefly B. tabaci, revealing the first global phylogeographic relationships for the group (Frohlich et al. 1999). The PCR reactions were performed in a 25 μl vol containing 1× Jumpstart REDTaq Ready-Mix (Sigma–Aldrich, Saint Louis, MI), primers (0.4 μM each), 20 ng of whitefly DNA and double distilled (dd) water in a 25 μl reaction volume. The PCR conditions for the 5′-end fragment were 35 cycles of 95°C for 60 s, 40°C for 60 s and 72°C for 90 s, with a final extension step of 7 min at 72°C. Cycling conditions for the 3′-end primers were accordingly, 95°C for 2 min, followed by 30 cycles of 95°C for 60 s, 52°C for 60 s and 72°C for 60 s, with a final extension for 5 min at 72°C. The expected size of the amplicons was confirmed by agarose gel (1%) electrophoresis in 1X TAE buffer, pH 8.0, containing 1X GelRed dye, at 100 V for 50 min. The PCR products of the expected sizes were cloned into the TA cloning plasmid vector pGEM-T Easy (Promega, Madison, WI), and sequenced bi-directionally using an ABI 3700 capillary sequencer at The University of Arizona Genomics Core facility (http://uagc.arl.arizona.edu/). Identification of B. tabaci Mitotypes From Pakistan The B. tabaci mitotypes were identified based on phylogenetic analysis of each of the COI fragments compared to a set of B. tabaci global reference sequences analogous to the respective COI fragment (Figs. 1 and 2). Reference sequences of both COI fragments representing the seven major phylogeographic clades of B. tabaci were downloaded from the NIH NCBI-GenBank database (Benson et al. 2012), or obtained from the laboratory COI DNA sequence database (JK Brown, unpublished data), containing sequences for adult B. tabaci field collections identified using the established taxonomic (morphological) characters (Takahashi 1936, Russell 1948) (courtesy, Mr. Ray Gill, California Department of Food and Agriculture, Sacramento, CA). Fig. 1. View largeDownload slide Phylogenetic reconstructions of the B. tabaci s.s.g. based on the (A) 5′- and (B) 3′-fragment. The samples used in this study were phylogenetically placed into three mitotypes within the Asia II phylogeographic clade. The nomenclature in parentheses corresponds to putative ‘species’ proposed by Dinsdale et al. (2010), or mitotypes, and bolded names are major phylogeographic clades, taken from Brown (2010). Only the Bayesian tree topologies are shown, and support values above branches correspond to bootstrap values / posterior probabilities at each node for the corresponding end, and the value placed below branches corresponds to bootstrap values / posterior probabilities at each node for the concatenated fragments. The greenhouse whitefly (Trialeurodes vaporariorum Westwood) was used as outgroup. Fig. 1. View largeDownload slide Phylogenetic reconstructions of the B. tabaci s.s.g. based on the (A) 5′- and (B) 3′-fragment. The samples used in this study were phylogenetically placed into three mitotypes within the Asia II phylogeographic clade. The nomenclature in parentheses corresponds to putative ‘species’ proposed by Dinsdale et al. (2010), or mitotypes, and bolded names are major phylogeographic clades, taken from Brown (2010). Only the Bayesian tree topologies are shown, and support values above branches correspond to bootstrap values / posterior probabilities at each node for the corresponding end, and the value placed below branches corresponds to bootstrap values / posterior probabilities at each node for the concatenated fragments. The greenhouse whitefly (Trialeurodes vaporariorum Westwood) was used as outgroup. Fig. 2. View largeDownload slide Observed mismatch distribution and fits to the expected model of demographic expansion. Mismatch distributions of individuals of the Asia II-1 haplotype (A) 5-fragment (raggedness statistic: 0.6242; r2: 0.0195) and (B) 3′-fragment (raggedness statistic: 0.6317; r2: 0.0068). Graphs C and D correspond to individuals of the Asia II-7 mitotype using the 5′- (raggedness statistic: 0.6848; r2: 0.0462) and 3′-fragment (raggedness statistic: 0.8955; r2: 0.0781), respectively. Fig. 2. View largeDownload slide Observed mismatch distribution and fits to the expected model of demographic expansion. Mismatch distributions of individuals of the Asia II-1 haplotype (A) 5-fragment (raggedness statistic: 0.6242; r2: 0.0195) and (B) 3′-fragment (raggedness statistic: 0.6317; r2: 0.0068). Graphs C and D correspond to individuals of the Asia II-7 mitotype using the 5′- (raggedness statistic: 0.6848; r2: 0.0462) and 3′-fragment (raggedness statistic: 0.8955; r2: 0.0781), respectively. Assembly and Alignment The COI sequences were assembled and manually edited using SeqMan Pro software available in DNASTAR Lasergene v8.0 (DNASTAR, Inc., Madison, WI). The assembled fragments were exported as FASTA files. The sequences were aligned using MUSCLE v3.8.31 (Edgar 2004) implemented in the Align Multiple Sequences tool, available in Mesquite v2.75 (Maddison and Maddison 2014), and the terminal gaps were treated as missing data. The sequence matrices consisting of three cloned amplicons per whitefly were used to detect nuclear insertions of the mitochondria (NUMTS) (see Screening for NUMTS). After removal of pseudogenes and chimeric sequences, a single mtCOI sequence per whitefly sample was selected for the identification of mitotypes and their spatial distribution pattern. To reconstruct the phylogenies, the identical or nearly identical haplotypes (99–100%) were identified using FABOX v1.41 (Villesen 2007), and redundant haplotype sequences were excluded from the alignment. Screening for NUMTS The aligned sequences were screened for NUMTS (or pseudogenes) according to the recommendations of Song et al. (2008). First, per above, cloning COI amplicons prior to DNA sequencing was considered a useful precaution against recovering NUMTS. Next, any sequence containing ambiguities, indels, stop codons, or singletons was not used in the analysis. This approach allowed multiple sequences cloned from each PCR product to be scrutinized, with those containing NUMTs standing out as divergent and/or containing stop codons. Then, based on the NJ tree reconstructed for the 3′- and 5′-end gene fragment, respectively, sequences that grouped in an unverifiable clade or formed a new clade (or sister clade) based on collection site and prior known histories of introductions, or if they grouped with an entirely different clade after several analyses using different groups of sequences, they were considered ‘suspect’ NUMTS and removed from the data set. Model Selection The optimal model of molecular evolution for the COI fragments was determined using jmodeltest v2.1.7 (Darriba et al. 2012), and a consensus of the Akaike Information Criterion (AIC), corrected AIC, and Bayesian Information Criterion (BIC). The HKY85, determined to be the best model of evolution (Hasegawa, 1985), was used for phylogenetic analysis of COI sequences, i.e., 3′-, 5′-, and the concatenated sequence consisting of both fragments for mitotype identification. Phylogenetic Analysis The mitochondria COI sequences were analyzed phylogenetically using Bayesian inference (BI) (MrBayes v3.2.5) (Huelsenbeck and Ronquist 2001, Ronquist et al. 2012), with four independent Markov Chain Monte Carlo (MCMC) runs consisting of four chains each, for 1 × 107 generations. Trees were sampled every 1,000 generations and the log-likelihood scores for the sampled points were plotted against the number of sampled generations using Tracer v1.6 (Rambaut et al. 2014) to verify stationarity of the chain and convergence. Effective sample size (ESS) values were tracked for each run until they exceeded at least 200 (Ronquist et al. 2012). The results were summarized using the commands, sump and sumt. The trees estimated before the chains reached a plateau were discarded as ‘burn-in’ trees, and represented the first 2 × 106 generations for each replicate. The 50% majority-rule consensus trees were drawn using FigTree v1.4.2 (http://tree.bio.ed.ac.uk/software/figtree/). The phylogenetic affiliation of each individual with a well-supported clade in the 3′-mtCOI tree was recorded and the geographic distribution of B. tabaci mitotypes was plotted using the ggplot2 in R package (Wickham 2009). Comparison Between the 3′-, 5′- and Concatenated COI Fragments The between-COI fragment sequence comparisons considered only sequences from samples for which an amplicon could be determined for both the 3′- and 5′-end of the gene, respectively. This had the effect of reducing the total number of whitefly samples used for this analysis from 829 to 538, alluding to the lower-than-expected specificity for both primer pairs, and suggesting greater than expected variability in the targeted regions. Sequence Polymorphisms The nt composition of the COI, number of segregating sites, and diversity were calculated using DNAsp v5.10.1 software (Librado and Rozas 2009). Nucleotide composition estimates do not provide sufficient information about substitution patterns, thus the ‘more informative’ ratio of transitions to transversions (Ti/Tv) was calculated using Mega v7.0 (Kumar et al. 2016). To determine whether the number and/or the type of substitutions may be biased by sequence composition differences, the degeneracy was assessed by calculating the ratio of 0-, 2- and 4-fold degenerate sites, using Mega v7.0 (Kumar et al. 2016). Intra- and Inter-Specific Sequence Divergence The intra- and inter-specific evolutionary pairwise distances were calculated for mitotype sequences analyzed by ML and BI, using the dist.dna function available in the APE package, implemented in R (Paradis et al. 2004). The HKY85 model of evolution (Hasegawa et al. 1985) was identified by jModel Test (Darriba et al. 2012) as the most optimal for all three alignment matrices, and used to correct distance estimations. Also, the Kimura 2-parameter (K2P) (Kimura 1980) was used to calculate intra- and inter-specific pairwise distances (Ross et al. 2008), with the ‘species delimitation’ plug-in available in Geneious (Masters et al. 2011). These analyses were carried out to determine whether the predicted mitotype composition in any given field population was expected to remain consistent, regardless of the different evolutionary estimates, or the COI gene fragment, e.g., 3′- or 5′-end, respectively. Finally, evolutionary estimates were used to compare the results based on previous predictions based on the K2P (Boykin et al. 2012) and HKY85 (Dinsdale et al. 2010) models, implemented for predicting species-cutoffs, at 1 and 3.5%, respectively, within the B. tabaci sibling species group. Phylogenetic Trees Reconstructed for the 3′- and 5′-COI Fragments To compare the phylogenies for the two COI regions and for the concatenated sequence consisting of both, sequences were aligned and subjected to Maximum Likelihood (ML) and BI analyses, using the HKY85 model of evolution, the best-fitting model for the different three alignments. The BI analysis was carried out, as described above. The ML tree was reconstructed using GARLI v2.1 (Zwickl 2006) with one thousand bootstrap iterations. The trees were rooted using the comparable COI fragments, respectively, from the greenhouse whitefly, Trialeurodes vaporariorum (Westwood) (GenBank Accession no. LN614547), as the outgroup. Consensus trees were drawn using FigTree v1.4.2 (http://tree.bio.ed.ac.uk/software/figtree/). The bootstrap values (bv) and the posterior probabilities (pp) were shown on the Bayesian tree reconstructed for each fragment separately, and for concatenated sequence consisting of both (Fig. 1). The ML tree topologies for the three alignments were statistically analyzed using the Shimodaira-Hasegawa test (SH test) (Shimodaira and Hasegawa 1999), available in PAUP v 4.0 (Swofford 2003). The SH test specifically compares the site log-likelihood differences for tree topologies. To test the hypothesis of similarity between trees, the null distribution of the test statistic was determined using nonparametric bootstrapping of the re-estimated log-likelihoods (RELL), and a two-tailed test, at a 95% confidence level. The three null hypotheses tested were: 1) there was no difference between 5′- and 3′-end tree topologies, 2) the 5′-fragment and concatenated COI sequences had the same topology, and 3) the topology of the 3′-COI fragment was the same as the concatenated tree. For each hypothesis, the SH test was performed using the DNA matrix for each pair of phylogenies, respectively. Population Genetics Analysis To determine whether a demographic shift could be detected based on the COI marker, the mutation frequency and the frequency of rare haplotypes within each population were evaluated using Dnasp v5.10.1 (Librado and Rozas 2009). The significance was tested by coalescence simulations using 10,000 replicates. The within-population mutation frequency of the 5′ and 3′-fragments and the 1269 base concatenated sequences, respectively, was calculated using Tajima’s D (Tajima 1989) and the R2 (Ramos-Onsins and Rozas 2002) test, respectively. The Tajima’s D compares the overall nt diversity in relation to diversity of segregating sites, while the R2 test calculates the differences for the number of singleton mutations and average number of nucleotides. The combined tests were selected to provide evidence for an ‘expansion hypothesis’, and in particular, the R2 was used due its suitability for analyzing small datasets (Ramirez-Soriano et al. 2008). The frequency of rare haplotypes within the B. tabaci populations was calculated using the Fu’s Fs statistic (Fu 1997), an approach that assumes no recombination at the locus of interest (Ramirez-Soriano et al. 2008). This test was carried out to provide a robust estimate of the number of ‘rare’ haplotypes, or those occurring in the population that represent recent mutations. Signals predicting change in population size were identified by plotting the frequencies of pairwise differences, referred to as ‘Mismatch distribution’, and the raggedness r statistic (Harpending 1994), where the extent of unimodality of the curve and low significant raggedness values are interpreted as evidence of population expansion. Finally, the most parsimonious network of haplotypes and their frequency were plotted using the Minimum Spanning Network (MSN) algorithm (Bandelt et al. 1999), implemented in PopART (http://popart.otago.ac.nz). The MSN analysis is expected to identify signals indicative of possible ‘demographic expansion.’ In an expansion scenario, the population structure assumes a ‘star-like shape’, representing predominant haplotypes surrounded less abundant ones, and those that differ by a few or rare mutations (Slatkin and Hudson 1991, Ruegg and Smith 2002, Ibáñez et al. 2011). Temporal Changes of DNA Polymorphisms in Populations of the Asia II-1 Mitotype To determine whether temporal patterns of demographic change observed for the Asia II-1 mitotype could be reconciled with the predicted genetic expansion followed by a recent bottleneck, the DNA polymorphisms were compared for the 3′-COI fragment obtained in this study (2012–2014) (n = 440), with the results reported for B. tabaci populations collected from the same locales in Pakistan during 1995–1998 (n = 6) (Simón et al. 2003) and 2005–2007 (n = 14) (Ahmed et al. 2010). Sequences corresponding to the 3′-COI end of the Asia II-1 mitotype were downloaded from GenBank using the batch entrez tool of the National Center for Biotechnology Information (NCBI) (https://www.ncbi.nlm.nih.gov/sites/batchentrez) and organized as comma-separated files (csv). The nt diversity (π), the diversity of haplotypes (Hd), and average number of nt differences (k) were calculated using DNAsp v6 (Rozas 2017). The polymorphism measures were compared statistically using a two-tailed t-test to evaluate the null hypothesis that the measures between two time periods were the same (Ho: µ1 = µ2). The following pairwise comparisons were analyzed: 1995–1998 by 2005–2007, 1995–1998 by 2012–2014, and 2005–2007 by 2012–2014. In addition, three measures of polymorphism per time period were compared using the χ2 goodness of fit test to assess the hypothesis that patterns of polymorphisms are consistent across the time periods studied. Finally, the rate of change across the three polymorphism indicators over time was calculated using a regression analysis. For all tests, significance was assessed at the 95% confidence level. Results Identification of B. tabaci Mitotypes From Pakistan The PCR amplification of the 5′- and 3′-COI fragments yielded a ~700 and ~850 bp amplicon, respectively, of which 618 and 651 bases were used in the analyses, after removing the primer sequences at both ends and carrying out quality trimming. The trimmed sequences were deposited in GenBank (Supp. Tables S1 and S2 for the 5- and 3-end, respectively). The global Bayesian analysis of the 5′- and 3′-fragment sequences showed that all of the haplotypes grouped in the major Asia II clade, with >95% pp support (Supp. Figs S1 and S2). Based on pairwise nt identity (Table 1) three subgroups (Dinsdale et al. 2010) were identified as Asia II-1, II-5, and II-7 (Supp. Figs. S1 and S2). The phylogenetic framework based on the 3′-COI fragment, resolved phylogenetic relationships that were consistent with previously published studies (Frohlich et al. 1999, Brown 2010, Dinsdale et al. 2010, Boykin et al. 2012) by grouping the genetic variants Asia II-1, II-2, II-3, and II-4 as a sister group to the group harboring mitotypes II-5 and -6, which themselves were sister groups to the II-7 mitotype. Also, the Asia II and Asia I major clades grouped monophyletically with reference mitotypes endemic to Australia, China (1, 2), and Italy, also consistent with previously defined monophyletic major clades (Brown 2010) (Supp. Fig. S2). In contrast, the relationships predicted by the 5′-COI sequence, were mostly in agreement with those based on the 3′-fragment; however, this analysis did not support the sistership between Asia II-1 and II-5 (Supp. Fig. S1), predicted by the 3′-fragment. Table 1. Sequence divergence and species delimitation measures for the three genetic variants found in Pakistan Haplotype Closest group Intra Dist (%) Intra HKY85 Inter Dist (%)a Inter HKY85b Clade supportc 5′-end Asia II-1 Asia II-7 0.215 0.578% 2.608 10.76% 1 Asia II-7 Asia II-1 0.116 0.074% 2.608 10.76% 1 Asia II-5 Asia II-7 0.081 0.065% 2.683 9.93% 1 Asia II-1 Asia II-5d - - - 11.65% - 3′-end Asia II-1 Asia II-7 0.142 0.674% 2.956 9.95% 1 Asia II-7 Asia II-1 0.068 0.016% 2.956 9.95% 1 Asia II-5 Asia II-1 0.239 0.375% 3.053 11.88% 1 Asia II-5 Asia II-7d - - - 12.50% - Concatenated Asia II-1 Asia II-7 0.17 0.627% 3.216 9.53% 1 Asia II-7 Asia II-1 0.071 0.044% 3.216 9.53% 1 Asia II-5 Asia II-7 0.147 0.200% 3.573 10.24% 1 Asia II-1 Asia II-5d - - - 10.71% - Haplotype Closest group Intra Dist (%) Intra HKY85 Inter Dist (%)a Inter HKY85b Clade supportc 5′-end Asia II-1 Asia II-7 0.215 0.578% 2.608 10.76% 1 Asia II-7 Asia II-1 0.116 0.074% 2.608 10.76% 1 Asia II-5 Asia II-7 0.081 0.065% 2.683 9.93% 1 Asia II-1 Asia II-5d - - - 11.65% - 3′-end Asia II-1 Asia II-7 0.142 0.674% 2.956 9.95% 1 Asia II-7 Asia II-1 0.068 0.016% 2.956 9.95% 1 Asia II-5 Asia II-1 0.239 0.375% 3.053 11.88% 1 Asia II-5 Asia II-7d - - - 12.50% - Concatenated Asia II-1 Asia II-7 0.17 0.627% 3.216 9.53% 1 Asia II-7 Asia II-1 0.071 0.044% 3.216 9.53% 1 Asia II-5 Asia II-7 0.147 0.200% 3.573 10.24% 1 Asia II-1 Asia II-5d - - - 10.71% - aInter-species tree pairwise distances calculated using the species delimitation determined using the Geneious plugin. bInter-species pairwise distances using the F84 model of evolution calculated in software available in the R package. cPosterior probability obtained from the Bayesian inference of the B. tabaci COI gene phylogeny. dGenetic variants are not sister clades and only HKY85 evolutionary inter-species distances was calculated for each group. View Large Table 1. Sequence divergence and species delimitation measures for the three genetic variants found in Pakistan Haplotype Closest group Intra Dist (%) Intra HKY85 Inter Dist (%)a Inter HKY85b Clade supportc 5′-end Asia II-1 Asia II-7 0.215 0.578% 2.608 10.76% 1 Asia II-7 Asia II-1 0.116 0.074% 2.608 10.76% 1 Asia II-5 Asia II-7 0.081 0.065% 2.683 9.93% 1 Asia II-1 Asia II-5d - - - 11.65% - 3′-end Asia II-1 Asia II-7 0.142 0.674% 2.956 9.95% 1 Asia II-7 Asia II-1 0.068 0.016% 2.956 9.95% 1 Asia II-5 Asia II-1 0.239 0.375% 3.053 11.88% 1 Asia II-5 Asia II-7d - - - 12.50% - Concatenated Asia II-1 Asia II-7 0.17 0.627% 3.216 9.53% 1 Asia II-7 Asia II-1 0.071 0.044% 3.216 9.53% 1 Asia II-5 Asia II-7 0.147 0.200% 3.573 10.24% 1 Asia II-1 Asia II-5d - - - 10.71% - Haplotype Closest group Intra Dist (%) Intra HKY85 Inter Dist (%)a Inter HKY85b Clade supportc 5′-end Asia II-1 Asia II-7 0.215 0.578% 2.608 10.76% 1 Asia II-7 Asia II-1 0.116 0.074% 2.608 10.76% 1 Asia II-5 Asia II-7 0.081 0.065% 2.683 9.93% 1 Asia II-1 Asia II-5d - - - 11.65% - 3′-end Asia II-1 Asia II-7 0.142 0.674% 2.956 9.95% 1 Asia II-7 Asia II-1 0.068 0.016% 2.956 9.95% 1 Asia II-5 Asia II-1 0.239 0.375% 3.053 11.88% 1 Asia II-5 Asia II-7d - - - 12.50% - Concatenated Asia II-1 Asia II-7 0.17 0.627% 3.216 9.53% 1 Asia II-7 Asia II-1 0.071 0.044% 3.216 9.53% 1 Asia II-5 Asia II-7 0.147 0.200% 3.573 10.24% 1 Asia II-1 Asia II-5d - - - 10.71% - aInter-species tree pairwise distances calculated using the species delimitation determined using the Geneious plugin. bInter-species pairwise distances using the F84 model of evolution calculated in software available in the R package. cPosterior probability obtained from the Bayesian inference of the B. tabaci COI gene phylogeny. dGenetic variants are not sister clades and only HKY85 evolutionary inter-species distances was calculated for each group. View Large The distribution and frequency of mitotypes varied depending on collection location. Based on the 3′-COI, approximately 83.8% (n = 588) of whiteflies were identified as mitotype Asia II-1, while, 15.7% grouped with Asia II-7 references (n = 110), and 0.5% clustered with Asia II-5 (n = 4). Similarly, based on the 5′-COI, 83.7% (n = 694) of whiteflies were identified as mitotype Asia II-1, while, 15.7% grouped with Asia II-7 references (n = 131), and 0.5% clustered with Asia II-5 (n = 4). The II-1 mitotype was the most abundant and broadly distributed, colonizing cotton and non-cotton hosts (Fig. 3), inferring it is a polyphagous phenotype, common among other Asia II mitotypes (Mound and Halsey 1978, Greathead 1986, Burban et al. 1992, Bedford et al. 1994, Qiu et al. 2007, Li et al. 2011). In addition to commercial cotton in the Punjab Province, Asia II-1 haplotypes were identified in the vicinity of Lahore, associated with herbaceous plant species grown in public and home gardens, small-scale vegetable plantings, and experimental research plot at the University of Punjab campus. In contrast, Asia II-7 and II-5 individuals were identified only in Lahore. Here, the II-7 mitotype colonized cotton, cucumber, okra, pepper, squash, and tomato plants, while the II-5 was associated with only cotton, cucumber, and squash plants, respectively, indicating that both were polyphagous while at the same time being more restricted in local geographic distribution, compared to the Asia II-1. Fig. 3. View largeDownload slide Map of Pakistan showing cotton-growing locales (source: USDA). Connectors indicate districts sampled and doughnut charts show the proportion of genetic group abundances in non-cotton (inner chart) and cotton (outer chart) hosts based on the 3′-COI fragment. Of the three locations studied, Lahore was the only locality where the three B. tabaci Asia II mitotype groups were found. Fig. 3. View largeDownload slide Map of Pakistan showing cotton-growing locales (source: USDA). Connectors indicate districts sampled and doughnut charts show the proportion of genetic group abundances in non-cotton (inner chart) and cotton (outer chart) hosts based on the 3′-COI fragment. Of the three locations studied, Lahore was the only locality where the three B. tabaci Asia II mitotype groups were found. Comparison of the 5′- and 3′- and Concatenated COI Fragments Phylogenetic Trees for the 5′- and 3′-COI Fragments The placements of the 5′- and 3′-COI fragments on the ML and BI trees, respectively, were consistent with the results of previous analyses placing Pakistan mitotypes in one of two distinct clades, Asia I and Asia II, two of seven recognized ‘major phylogeographical clades’, with pp and bv of >99% and >70%, respectively (Fig. 1). However, the predicted sister-relationships were either in conflict and/or unresolved regardless of the method used for phylogenetic analysis. Within the Asia II major clade, analysis of the 5′-fragment provided support for a sister relationship between the II-1 and II-7 Pakistan mitotypes by the BI but not ML. This result was consistent with the lack of robust support for the concatenated sequence analysis, with a posterior probability of 77% (BI) and bv of <50% (Fig. 1). By comparison the BI analysis of the 3′-end COI sequence predicted a sister relationship between II-1 and II-5; however, it was only weakly supported at 66% pp. The B and Q mitotypes were predicted to be sister groups and relationships among these members of the North Africa-Mediterranean-Middle Eastern (NAF_MED_ME) major clade were well-supported, irrespective of COI fragment. Of the three BI analyses, only the 3′-COI and the concatenated sequences predicted well-supported sister association between the NAF_MED_ME and Asia II clades, at 99% pp, compared to the 5′-COI sequence tree (Fig. 1). Finally, both ML and BI analyses supported the two ‘American Tropics’ (New World) monophylies on the 5′-, 3′-, and concatenated trees (Fig. 1). Interestingly, their sistership to the Asia II and North Africa-Mediterranean-Middle East clades, respectively, was supported by BI analysis of the 3′-, 5′-, and concatenated sequences, whereas, based on ML analysis only the concatenated fragment supported this association, at 75% bootstrap (Fig. 1). For the major Sub-Saharan Africa (SSA) clade, three sister mitotypes (SSA1-SSA3) were resolved based on BI and ML analysis of both fragments and the concatenated sequence, each of which were well-supported (pp > 95%, bv > 75). However, the monophyly of the major clade (SSA) that contains the three mitotypes, was well-supported by BI analysis of both fragments, and only weakly supported by the ML inference for the 5′- and 3′-fragments, at 63 and 60 bv, respectively. Within the major clade (SSA), a discrepancy was observed in the BI analysis, which predicted a sister relationship between the SSA-1 and SSA-2 mitotypes, with robust support for the 3′-, but not for the 5′-end or the concatenated sequence (Fig. 1). Results of the SH tree topology test indicated that the ‘best trees’ inferred by ML for the 5′- and 3′- fragments were not significantly different (P = 0.175) (Table 2). However, a significant difference was observed in the trees reconstructed for the 5′- and 3′-end-, and concatenated sequence, respectively (P < 0.05) (Table 2). Table 2. Shimodaira-Hasegawa test using RELL bootstrap for comparison of phylogenetic tree topologies between COI fragments and the concatenated sequence Diff -ln L P-value 5′-end vs 3′-end 349.3876 0.175 5′-end vs concatenated 1139.633 0.001* 3′-end vs concatenated 790.2452 0.018* Diff -ln L P-value 5′-end vs 3′-end 349.3876 0.175 5′-end vs concatenated 1139.633 0.001* 3′-end vs concatenated 790.2452 0.018* *P < 0.05, differences in topologies are significant at the 0.05 level. View Large Table 2. Shimodaira-Hasegawa test using RELL bootstrap for comparison of phylogenetic tree topologies between COI fragments and the concatenated sequence Diff -ln L P-value 5′-end vs 3′-end 349.3876 0.175 5′-end vs concatenated 1139.633 0.001* 3′-end vs concatenated 790.2452 0.018* Diff -ln L P-value 5′-end vs 3′-end 349.3876 0.175 5′-end vs concatenated 1139.633 0.001* 3′-end vs concatenated 790.2452 0.018* *P < 0.05, differences in topologies are significant at the 0.05 level. View Large The unexpected placement of an unverified Asian mitotype from Pakistan, GenBank Accession no. KJ709461 (Ashfaq 2014), with the ‘presumed monophagous’ sweet potato mitotype from Uganda (UG-SwPot) (Accession no. AY057174) was supported by the ML concatenated tree, but not by the 5′- or 3′-sequence trees. However, the BI trees reconstructed for the separate COI fragments and concatenated version showed robust support (Fig. 1). Even so, in previous studies the UG-SwPot mitotype has grouped separately within the sub-Saharan African clade, and therefore is a distinct mitotype uniquely associated with sweet potato in east Africa (Legg et al. 2002, Boykin et al. 2007, see references in Brown et al. 2010, Gill and Brown 2010, Firdaus et al. 2013, Legg et al. 2014). Given that a sequence is available for only a single individual of this uncharacterized whitefly mitotype and that the origin of its closest relative is Africa instead of Asia, as might be expected, the possibility that the sequence is a NUMTs (Song et al. 2008, Tay et al. 2017) requires further investigation. Nucleotide Sequence Polymorphisms The proportions of T, C, A, and G bases associated with the 5′-COI fragment were 44.3%, 14.2%, 23.2%, and 18.3% respectively, compared to the 3′-end sequence, at 43.8%, 13.1%, 22.9%, 20.3%, respectively. The percentage GC content of the two COI fragments was nearly identical, at 32% and 33% respectively, with differences occurring in the third codon position (Table 3). Of the two fragments, the 3′-sequence harbored more segregating sites and mutations, while the 5′-fragment exhibited greater overall greater nt diversity (π). The nt diversity (π) in the second codon position of the 5′-end COI was approximately three times greater compared to the 3′-end sequence (Table 4). And, the Ti/Tv ratio for the 5′- compared to the 3′-sequence was nearly twice as great, at 5.9 and 2.93, respectively (Table 5). Table 3. Nucleotide composition of the 3′- and 5′-end COI fragments Codon position 5′-end 3′-end T (%) A (%) G (%) C (%) T (%) A (%) G (%) C (%) Overall 44 23 18 14 44 23 20 13  First 31 29 25 15 30 29 26 15  Second 50 13 17 20 47 16 20 17  Third 52 27 13 7 54 24 15 7 Codon position 5′-end 3′-end T (%) A (%) G (%) C (%) T (%) A (%) G (%) C (%) Overall 44 23 18 14 44 23 20 13  First 31 29 25 15 30 29 26 15  Second 50 13 17 20 47 16 20 17  Third 52 27 13 7 54 24 15 7 View Large Table 3. Nucleotide composition of the 3′- and 5′-end COI fragments Codon position 5′-end 3′-end T (%) A (%) G (%) C (%) T (%) A (%) G (%) C (%) Overall 44 23 18 14 44 23 20 13  First 31 29 25 15 30 29 26 15  Second 50 13 17 20 47 16 20 17  Third 52 27 13 7 54 24 15 7 Codon position 5′-end 3′-end T (%) A (%) G (%) C (%) T (%) A (%) G (%) C (%) Overall 44 23 18 14 44 23 20 13  First 31 29 25 15 30 29 26 15  Second 50 13 17 20 47 16 20 17  Third 52 27 13 7 54 24 15 7 View Large Table 4. The DNA polymorphisms and nucleotide diversity for the 5′- and 3′-COI fragments Parameter Codon position 5′-end (n = 538) 3′-end (n = 538) S Overall 71 98  First 13 15  Second 11 12  Third 47 71 Hd Overall 0.398 ± 0.001 0.420 ± 0.0011**  First 0.340 ± 0.001 0.349 ± 0.001**  Second 0.327 ± 0.001 0.358 ± 0.001**  Third 0.380 ± 0.001** 0.371 ± 0.0011 π Overall 0.028 ± 0.0001** 0.026 ± 0.0001  First 0.016 ± 0.0001** 0.006 ± 0.00002  Second 0.0027 ± 9.6E-06** 0.002 ± 0.00001  Third 0.068 ± 0.0002 0.072 ± 0.0002** θ Overall 0.027 ± 0.0003 0.028 ± 0.00025*  First 0.0138 ± 0.0002** 0.012 ± 0.0002  Second 0.0122 ± 0.0002** 0.010 ± 0.00015  Third 0.057 ± 0.001 0.063 ± 0.0006** κ Overall 10.784 13.264  First 2.266 1.034  Second 0.352 0.395  Third 8.166 13.264 Parameter Codon position 5′-end (n = 538) 3′-end (n = 538) S Overall 71 98  First 13 15  Second 11 12  Third 47 71 Hd Overall 0.398 ± 0.001 0.420 ± 0.0011**  First 0.340 ± 0.001 0.349 ± 0.001**  Second 0.327 ± 0.001 0.358 ± 0.001**  Third 0.380 ± 0.001** 0.371 ± 0.0011 π Overall 0.028 ± 0.0001** 0.026 ± 0.0001  First 0.016 ± 0.0001** 0.006 ± 0.00002  Second 0.0027 ± 9.6E-06** 0.002 ± 0.00001  Third 0.068 ± 0.0002 0.072 ± 0.0002** θ Overall 0.027 ± 0.0003 0.028 ± 0.00025*  First 0.0138 ± 0.0002** 0.012 ± 0.0002  Second 0.0122 ± 0.0002** 0.010 ± 0.00015  Third 0.057 ± 0.001 0.063 ± 0.0006** κ Overall 10.784 13.264  First 2.266 1.034  Second 0.352 0.395  Third 8.166 13.264 S, Number of segregating sites; Hd, Haplotype diversity; π, Nucleotide diversity; θ, Nucleotide diversity from S; κ, Average number of nucleotide differences. Diversity values correspond to the mean ± the SEM. *The mean difference is significant at the 0.01 level; **The mean difference is significant at the 0.001 level. Asterisks stand by the highest mean. View Large Table 4. The DNA polymorphisms and nucleotide diversity for the 5′- and 3′-COI fragments Parameter Codon position 5′-end (n = 538) 3′-end (n = 538) S Overall 71 98  First 13 15  Second 11 12  Third 47 71 Hd Overall 0.398 ± 0.001 0.420 ± 0.0011**  First 0.340 ± 0.001 0.349 ± 0.001**  Second 0.327 ± 0.001 0.358 ± 0.001**  Third 0.380 ± 0.001** 0.371 ± 0.0011 π Overall 0.028 ± 0.0001** 0.026 ± 0.0001  First 0.016 ± 0.0001** 0.006 ± 0.00002  Second 0.0027 ± 9.6E-06** 0.002 ± 0.00001  Third 0.068 ± 0.0002 0.072 ± 0.0002** θ Overall 0.027 ± 0.0003 0.028 ± 0.00025*  First 0.0138 ± 0.0002** 0.012 ± 0.0002  Second 0.0122 ± 0.0002** 0.010 ± 0.00015  Third 0.057 ± 0.001 0.063 ± 0.0006** κ Overall 10.784 13.264  First 2.266 1.034  Second 0.352 0.395  Third 8.166 13.264 Parameter Codon position 5′-end (n = 538) 3′-end (n = 538) S Overall 71 98  First 13 15  Second 11 12  Third 47 71 Hd Overall 0.398 ± 0.001 0.420 ± 0.0011**  First 0.340 ± 0.001 0.349 ± 0.001**  Second 0.327 ± 0.001 0.358 ± 0.001**  Third 0.380 ± 0.001** 0.371 ± 0.0011 π Overall 0.028 ± 0.0001** 0.026 ± 0.0001  First 0.016 ± 0.0001** 0.006 ± 0.00002  Second 0.0027 ± 9.6E-06** 0.002 ± 0.00001  Third 0.068 ± 0.0002 0.072 ± 0.0002** θ Overall 0.027 ± 0.0003 0.028 ± 0.00025*  First 0.0138 ± 0.0002** 0.012 ± 0.0002  Second 0.0122 ± 0.0002** 0.010 ± 0.00015  Third 0.057 ± 0.001 0.063 ± 0.0006** κ Overall 10.784 13.264  First 2.266 1.034  Second 0.352 0.395  Third 8.166 13.264 S, Number of segregating sites; Hd, Haplotype diversity; π, Nucleotide diversity; θ, Nucleotide diversity from S; κ, Average number of nucleotide differences. Diversity values correspond to the mean ± the SEM. *The mean difference is significant at the 0.01 level; **The mean difference is significant at the 0.001 level. Asterisks stand by the highest mean. View Large Intra- and Inter-Specific Sequence Divergence The average intraspecific HKY85 pairwise divergence for the 5′-COI-end sequence among the Asia II subgroups, was 0.578% for Asia II-1 (n = 440), 0.074% for II-7 (n = 94), and 0.065% for II-5 (n = 4). In contrast, pairwise divergence for the 3′-sequences was 0.674% for Asia II-1 (n = 440), and 0.016% and 0.375% for II-7 (n = 94) and II-5 (n = 4), respectively (Table 1). Table 5. Evolutionary parameters for the 3′- and 5′-COI fragments estimated from 538 sequences and three mitotypes Codon position 5′-end 3′-end Ti/Tv ratio Overall 5.9 2.93  First 4.81 1.85  Second 0.99 1.32  Third 5.41 4.22 0-fold degenerated sites ratio 0.62 0.63 Twofold degenerated sites ratio 0.17 0.18 Fourfold degenerated sites ratio 0.15 0.14 Codon position 5′-end 3′-end Ti/Tv ratio Overall 5.9 2.93  First 4.81 1.85  Second 0.99 1.32  Third 5.41 4.22 0-fold degenerated sites ratio 0.62 0.63 Twofold degenerated sites ratio 0.17 0.18 Fourfold degenerated sites ratio 0.15 0.14 View Large Table 5. Evolutionary parameters for the 3′- and 5′-COI fragments estimated from 538 sequences and three mitotypes Codon position 5′-end 3′-end Ti/Tv ratio Overall 5.9 2.93  First 4.81 1.85  Second 0.99 1.32  Third 5.41 4.22 0-fold degenerated sites ratio 0.62 0.63 Twofold degenerated sites ratio 0.17 0.18 Fourfold degenerated sites ratio 0.15 0.14 Codon position 5′-end 3′-end Ti/Tv ratio Overall 5.9 2.93  First 4.81 1.85  Second 0.99 1.32  Third 5.41 4.22 0-fold degenerated sites ratio 0.62 0.63 Twofold degenerated sites ratio 0.17 0.18 Fourfold degenerated sites ratio 0.15 0.14 View Large The evolutionary divergences for the 3′- and 5′-sequences differed among the whitefly populations. Using the HKY85 model of evolution, the intra-species divergence for the 5′-sequence of the II-7 group was five times greater than for the 3′-sequence. By comparison, for the II-5 mitotype 3′-end, the genetic divergence was twice as great as for the 5′-end sequence. Similar patterns were observed for the K2P model (Ross et al. 2008) in which the 5′-fragment diverged by 0.215% for II-1 (n = 440), 0.116% for Asia II-7 (n = 94), and 0.081% for II-5 (n = 4). By comparison, the 3′-fragment for II-1 diverged by 0.142% (n = 440), whereas, II-7 and II-5 diverged by 0.068% (n = 94) and 0.239% (n = 4), respectively (Table 1). Also, the analysis of the concatenated sequences (3′- and 5′-ends) using the HKY85 model indicated a range of divergence from 0.044 to 0.627%, compared to 0.071 to 0.17% based on the K2P model. Finally, the inter-species HKY85 distances for the 5′-end sequences ranged from 9.93% to 11.65%, while the 3′-end showed slightly higher divergences ranging between 9.95 and 12.50%. Population Genetics of Mitotypes in Pakistan The distribution of frequencies for pairwise comparisons of the 3′- and 5′-COI fragments, individually, and concatenated are summarized in Fig. 2. Based on these results the Asia II-1 and II-7 populations had unimodal mismatch distributions, a distribution that is strongly associated with genetic expansion at the population level. Because of the small sample size it was not possible to determine the mismatch distribution and population genetics measures for Asia II-5 (n = 4). However, the smoothness of the curves, based on the raggedness r statistic, indicated no significant differences among the populations analyzed here. The Tajima’s D values for Asia II-1 and Asia II-7 were significantly negative for the 5′-or the 3′-ends, and the concatenated sequence, respectively. The negative Tajima’s D value reflected low variation, compared to the expected level, leading to the hypothesis that this population underwent a demographic expansion and/or ‘recovery’, commonly observed following bottleneck events, regardless of mechanism (Table 6). In further support, the Ramos-Onsin and Rozas’s R2 values, which were at near-zero, supported a population growth hypothesis for both Asia II-1 and Asia II-7 (Table 6). Also, Fu’s Fs statistic, which was significantly negative, supported the ‘population growth’ hypothesis (Table 6). The statistical values and levels of significance were about the same for the 5′- and 3′-COI sequences, or the concatenated sequence, respectively (Table 6). Table 6. Statistical tests assessing population expansion of Asia II-1 and Asia II-7 mitotypes based on the 3′-, 5′-, and concatenated fragments 5′-end 3′-end Concatenated R2 Asia_II-1 0.0122* 0.0171* 0.0124* Asia II-7 0.0485* 0.0844* 0.0538* Tajima’s D Asia II-1 −2.28947** −2.61804*** −2.69952*** Asia II-7 −2.49512** −1.91738* −2.60292*** Fu’s Fs Asia II-1 −33.603* −44.313* −83.631* Asia II-7 −6.948* −2.243* −7.156* 5′-end 3′-end Concatenated R2 Asia_II-1 0.0122* 0.0171* 0.0124* Asia II-7 0.0485* 0.0844* 0.0538* Tajima’s D Asia II-1 −2.28947** −2.61804*** −2.69952*** Asia II-7 −2.49512** −1.91738* −2.60292*** Fu’s Fs Asia II-1 −33.603* −44.313* −83.631* Asia II-7 −6.948* −2.243* −7.156* *The mean difference is significant at the 0.05 level; **The mean difference is significant at the 0.01 level; ***The mean difference is significant at the 0.001 level. View Large Table 6. Statistical tests assessing population expansion of Asia II-1 and Asia II-7 mitotypes based on the 3′-, 5′-, and concatenated fragments 5′-end 3′-end Concatenated R2 Asia_II-1 0.0122* 0.0171* 0.0124* Asia II-7 0.0485* 0.0844* 0.0538* Tajima’s D Asia II-1 −2.28947** −2.61804*** −2.69952*** Asia II-7 −2.49512** −1.91738* −2.60292*** Fu’s Fs Asia II-1 −33.603* −44.313* −83.631* Asia II-7 −6.948* −2.243* −7.156* 5′-end 3′-end Concatenated R2 Asia_II-1 0.0122* 0.0171* 0.0124* Asia II-7 0.0485* 0.0844* 0.0538* Tajima’s D Asia II-1 −2.28947** −2.61804*** −2.69952*** Asia II-7 −2.49512** −1.91738* −2.60292*** Fu’s Fs Asia II-1 −33.603* −44.313* −83.631* Asia II-7 −6.948* −2.243* −7.156* *The mean difference is significant at the 0.05 level; **The mean difference is significant at the 0.01 level; ***The mean difference is significant at the 0.001 level. View Large The results of the MSN analysis indicated a likely demographic expansion, as reflected by the star-like structures that depict one predominant haplotype surrounded by several less abundant or ‘rare’ haplotypes harboring a few mutations, low-level sequence divergence and a high frequency of unique mutations (Mirol et al. 2008). This result agreed with the low, observed variation predicted by Tajima’s D test (Table 6). For the 5′-end COI sequence the predominant network haplotypes represented 94% of Asia II-1 individuals, and 96% of Asia II-7 sequences (Fig. 4A). However, for the 3′-COI marker, 93% of Asia II-1 (1% less than in the 5′-end) and 98% of Asia II-7 (2% more than in the 5′-end) accounted for the common haplotypes (Fig. 4B). These results were consistent with number of segregating sites for Asia II-1 and II-7, based on the nt diversity and sequence divergence measures, indicating more variation in the 3′- than 5′-sequence. These results collectively confirmed that the 3′- and 5′-COI gene fragments were not interchangeable, and that when taken together, e.g., concatenated sequence (Fig. 4C), they contributed more valuable information, or were more informative together than when analyzed separately. Fig. 4. View largeDownload slide Minimum Spanning Network (ɛ = 0) based on 618 bp of the 5′ end (A), 651 bp of the 3′-end (B), and 1227 bp of the concatenated 5′- and 3′-fragments (C) of the COI gene. Circles represent haplotypes, and the size is proportional to frequency. The star-like structures showing a predominant haplotype, surrounded by rare haplotypes that have accumulated 1-5 mutations, associated with Asia II-1 (left) and-7 (right), in support of demographic expansion. Small vertical bars across connectors show the number of mutations/steps between haplotypes. One hypothetical most parsimonious route is shown between mitotypes Asia II-1 and −7. Pie charts show the proportion of plant hosts, as cotton or non-cotton for each individual sample. Fig. 4. View largeDownload slide Minimum Spanning Network (ɛ = 0) based on 618 bp of the 5′ end (A), 651 bp of the 3′-end (B), and 1227 bp of the concatenated 5′- and 3′-fragments (C) of the COI gene. Circles represent haplotypes, and the size is proportional to frequency. The star-like structures showing a predominant haplotype, surrounded by rare haplotypes that have accumulated 1-5 mutations, associated with Asia II-1 (left) and-7 (right), in support of demographic expansion. Small vertical bars across connectors show the number of mutations/steps between haplotypes. One hypothetical most parsimonious route is shown between mitotypes Asia II-1 and −7. Pie charts show the proportion of plant hosts, as cotton or non-cotton for each individual sample. Temporal Changes of DNA Polymorphism in Populations of the Asia II-1 Mitotype Changes in DNA polymorphisms in the mtCOI for B. tabaci sampled from the Pakistani cotton over the 1995–2014 timeframe (Fig. 5) were observed in that haplotype diversity and average number of nt differences remained unchanged during the earlier sampling periods but dropped significantly for 2012–2014 population. In contrast, changes in nt diversity began to show significant differences earlier in the 2005–2007 population and continued to decline in the subsequent sampling period. A significant, negative correlation between the three polymorphic indicators was statistically supported by the regression analysis, in that, a significant progressive reduction in diversity occurred over a short time. The rate of change inferred from the r2 value for the haplotype diversity index showed a less dramatic deviation than for the other two indicators (r2 > 0.9), likely attributable to the nature of the calculation of this index that compares sequences as a whole and not site by site. The changes in the average number of nt differences (k) and the nt diversity showed strong r2 values, supporting the dramatic change in diversity for the last sampling period. Additionally, the χ2 goodness of fit test allowed the rejection of the null hypothesis (P < 0.05) that the pattern of polymorphism indicators remained unchanged over the three annual sampling periods. These results were consistent with a ‘recovery’ scenario, in that the population was observed to undergo genetic expansion following an apparent recent bottleneck. Such a pattern has often been associated with overuse of pesticides, resulting in the development of insect pest resistance, with several well-known examples documented in whitefly heavily exposed to organophosphate chemistries (Cahill et al. 1995; Ahmad et al. 2002, 2010). Fig. 5. View largeDownload slide Temporal changes in DNA polymorphism of Asia II-1 populations sampled from the cotton-growing area in Pakistan. Different letters represent statistical differences at the 95% confidence level (P < 0.05) based on a two-tailed t-test. The rate of change is represented by the regression line and the significant r2 values. Fig. 5. View largeDownload slide Temporal changes in DNA polymorphism of Asia II-1 populations sampled from the cotton-growing area in Pakistan. Different letters represent statistical differences at the 95% confidence level (P < 0.05) based on a two-tailed t-test. The rate of change is represented by the regression line and the significant r2 values. Discussion The cotton leaf curl disease (CLCuD) has resulted in huge economic losses in cotton production in Pakistan (Farooq et al. 2014), however, shifts in whitefly vector dynamics and their association with the spread of CLCuMuV and later the resistance-breaking CLCuKoV-Bur strain in the Punjab and Sindh provinces (Rajagopalan et al. 2012, see references in Ashfaq et al. 2014, Sattar et al. 2017) has not been investigated. In this study, the abundance and distribution of whitefly haplotypes, niche association, and changes in the structure of haplotype diversity were studied for the predominant B. tabaci mitotypes in cotton plantings in Pakistan. In particular, based on population genetics analyses, the Asia II-1 mitotype was shown to have ‘recovered’ following a recent bottleneck, an event that has dramatically reduced the genetic variability within Asia II-1, and the diversity if B. tabaci mitotypes in cotton crops in the locales studied here. The patterns of DNA polymorphisms associated with Asia II-1 collections collected in three time periods since 1990, somewhat suggested that selection of this mitotype was associated with pesticide use. B. tabaci Mitotype Diversity, Spatial Distribution, and Host Association With respect to mitotype characterization and spatial distribution patterns, this study, together with previous reports, provides evidence that B. tabaci in Pakistan are represented by genetic and phenotypic variants that occupy different geographical and host range niches. Recently, Masood et al. (2017) reported that Asia II-7 was restricted to the northeast of Punjab, where it occurred at a very low frequency, whereas, Asia II-5 was sporadically distributed in cotton in the northern and southern Punjab. In this study, three B. tabaci mitotypes were identified, Asia II-1, −5, and −7 (Dinsdale et al. 2010, Boykin et al. 2012), associated with cotton and/or non-cotton plant host species in Pakistan. In 1992, Bedford et al. (1994) reported the biological, genetic, and morphological characteristics for the whitefly colony, PC92 that grouped within the Asia II-1, established from collections from Multan cotton fields (Bedford et al. 1994, Rosell et al. 1997, Brown 2000). At about the same time, members of the Asia II-1 were found colonizing watermelon plants in northern India (IW) in 1991, and Nepal (NeW) in 1994, from collection sites at approximately the same latitude. The presence of these mitotype group in Sindh, a southern province of Pakistan where the group was not known to colonize (Ahmed et al. 2010, 2011), was recently documented in 2013 (Ashfaq et al. 2014), suggesting that Asia II-1 have probably been established in the region for some time, before their southerly expansion that probably began during the year 2000, or slightly later. These shifts in whitefly distribution in Pakistan are time-wise consistent with the initial CLCuMV outbreak and virus spread throughout the Punjab province, and the emergence and spread of the recombinant CLCuKoV-Bu strain, which began in ~2002 and continues to the present (Amrao et al. 2010). In addition to the haplotypes that group in clades comprising the different mitotypes groups endemic to Asia, the B mitotype has been documented in Pakistan. The earliest reports of B mitotype (North Africa-Mediterranean-Middle East major clade) introductions to the Indian subcontinent were in India in 1999 (Banks et al. 2001) and Pakistan in 1998 Simón et al. (2003), respectively. During 2011–2012, two independent studies confirmed that the B mitotype had become established in Pakistan (Ahmed et al. 2011, Hameed et al. 2012), and indeed, subsequent field surveys confirmed endemic Asia I and Asia II clades, as well as the exotic B mitotype (Ashfaq et al. 2014, Masood et al. 2017). However, in Pakistan, the distribution of the B type has been reported to be limited to the southern province of Sindh (Ahmed et al. 2010, 2011; Ashfaq et al. 2014; Masood et al. 2017), suggesting that haplotypes from this group was not involved in virus spread during the outbreaks occurring in Multan and Burewala in the Punjab province (Mansoor 2003a), a location where the Asia II-1 is the predominant type. The extent of polyphagy or oligophagy, namely, the particular host range, as well as host preference, among the different mitotypes of the B. tabaci complex (Mound 1978, Cock 1993, Mound 1963, Bedford et al. 1994, Fishpool and Burban 1994) are factors that can be considered in the progressive expansion hypothesis observed here for the putative Asia II-1. Although not essential, polyphagy among expanding populations would be expected to enhance fitness. This is consistent with the broad host range and wide distribution known for Asia II-1 collected in China and India during the 1990s to mid 2000s (Qiu et al. 2007, 2011) that reported this mitotype group in cassava, cotton, lantana, pumpkin, sunflower, and tobacco. Similarly, in this study, members of the Asia II-1 were collected in cotton, cucumber, okra, pepper, squash, and tomato. In addition to the polyphagous behavior of the Asia II-1, the mitotype group Asia ll-7 has also been reported to be polyphagous in that it has been reported from Hibiscus, laurel, poinsettia, and tomato in China (Qiu et al. 2011), and herein, from cotton, cucumber, okra, pepper, squash, and tomato plants. These observations, suggest a role for polyphagy in the recent population dispersal to the south, albeit, this trait is only one of several that could be involved, together with insecticide resistance, factors related to the chemical ecology of virus-infected plants and/or virus-acquisition that are expected to alter B. tabaci behavior and potentially fitness, among others (Byrne and Bellows 1991, Costa et al. 1991, Liu et al 2007, Singer 2008). Population Genetics of Asia II-1 and II-7 Populations Here, several lines of evidence are provided indicating that both genetic and demographic expansion have recently occurred in the II-1 and II-7 populations, a conclusion based on Tajima’s D, Fu’s Fs, and Ramos-Onsins and Rozas’s R2 tests. The haplotype structure observed can be explained by an excess of low-frequency polymorphisms and the unimodality of pairwise distribution associated with the II-1 and II-7 mitotypes (Rogers and Harpending 1992, Harpending 1994). Although the observed low sequence diversity has been hypothesized to result from selective sweeps associated with certain whitefly-endosymbiotic bacteria (Hurst and Jiggins 2005), a demographic expansion/recovery followed by a recent bottleneck may likely explain the apparent (predicted) rapid re-distribution of the Asia II-1 mitotype in the Punjab and Sindh provinces in Pakistan. This hypothesis is supported by the dramatic change of DNA polymorphisms observed over the last two decades, whereby the selection pressure exerted by the overuse of pesticides provoked a loss in intraspecific diversity and gain of resistance to widely used chemical groups in the whitefly pest and virus vector (Cahill et al. 1995) that aided in an outcome in which certain haplotypes of the Asia II-1 thrived over susceptible individuals, and specifically over haplotypes of Asia II-7 that exhibited higher susceptibility to organophosphates than other Asian mitotypes (Naveen et al. 2017). Although the putative Asia II-7 and II-1 mitotypes appear to have been exposed to similar insecticides in the cotton-growing locales, that the Asia II-1 type thrived in the cotton-growing locales, while II-7 did not, may explain the persistence of the II-7 near Lahore where pesticide use is overall less or non-existent in urban and small vegetable plantings. The temporal changes in patterns in polymorphism are consistent with the pesticide use patterns documented in Pakistan. Following the shift from older cotton varieties in Pakistan to the S12 cotton genotype, a hairy (hirsute) whitefly-susceptible variety, whitefly populations increased substantially for the first time in the Punjab province, leading to the first cotton leaf curl disease outbreak that began in approximately 1990 (Briddon and Markham 2000). This led to the heavy use of organophosphates, and to the rapid development of high levels of resistance in the endemic Asia II B. tabaci, possibly, Asia II-1. Pressure exerted by the overuse of organophosphates was alleviated during the next decade when a number of new, more selective insecticides were introduced, together with staggered use to manage resistance, including insect-growth regulators, neonicotinoids, and pyrethroids, among others (Ahmad et al. 2002). One example of such a compound was imidacloprid (neonicotinoid) to which minimal or no resistance was documented during 1996 to 2011 among B. tabaci populations (Ahmad and Khan 2017). Accordingly, it appears that as the result of these whitefly management practices, little or no discernable change occurred in haplotypic structure of Asia II-1 from the mid-1990s to the early 2000s. However, by 2007, following the outbreak of the recombinant CLCuKoV-Bu (Amrao et al. 2010), hard pesticide chemistry-use increased substantially in cotton to combat virus spread, by ~12% (Tariq et al. 2007), and resistance continued to rise (Ahmad et al. 2010, Ahmad and Khan 2017). The consequence of the change in pesticide use patterns appears to have contributed importantly to the loss in haplotype diversity that by the 2012–2014 period, declined dramatically (Fig. 5). Even though the apparent reduction in diversity may have been over- or under-estimated due to the relatively small sample size, the well-documented shifts in agronomic practices supported a trend toward a temporal reduction in diversity that preceded the ‘recovery’ predicted in this study. Other phenomena such as hormoligosis, together with destabilization of a pest population away from equilibrium resulting by the concomitant decline in beneficial insects due to pesticide use (Dittrich 1987), appears to also have contributed importantly to the observed Asia II-1 population expansion. Based on historical and current documentation, it can be hypothesized that the expansion of the Asia II-1 population has contributed twice to begomovirus outbreaks of in cotton-vegetable cropping systems, first in ~1990, followed by emergence and spread of CLCuMuV and other begomoviral species/strains (Brown 2017), and again during 2002-onward, with the release of CLCuMuV-resistant germplasm, which fostered the CLCuKoV-Bu outbreak that has reached pandemic proportions throughout most Pakistan cotton-growing regions and northern India (Punjab). This emergence of the recombinant pandemic CLCuKoV-Bu strain appears to have emerged as a direct result of the widespread cultivation of a CLCuMuV-resistant cotton variety developed to combat the 1990 outbreak (Briddon et al. 2000, Mansoor et al. 2003a). A pattern of initial disease outbreak caused by a complex of strains and species, followed by the emergence of a recombinant begomovirus, driven in part by a population expansion of the whitefly vector, is reminiscent of the phylodynamics surrounding the cassava mosaic disease (CMD) outbreak, which also involved the expansion of a whitefly B. tabaci vector mitotype in Sub-Saharan Africa (Legg et al. 2002, 2014). Although previous studies have considered the association between the increased B. tabaci population size in cotton-vegetable cropping systems in Pakistan and cotton leaf curl disease outbreaks, and have identified several B. tabaci mitotypes from 1990 onward (Costa et al. 1993; Bedford et al. 1994; Frohlich et al. 1999; Simón 2003; Ahmed et al. 2010, 2011; Masood et al. 2017), this report offers the first timetable of and evidence for a genetic expansion hypothesis involving the Asia II-1 mitotype. It is well known that factors such as introductions or pesticide-induced upsurges of exotic plant viruses and/or whitefly vectors can alter the phylodynamics of plant viruses, particularly in cropping systems (Brown and Bird 1992, Brown 2007). In some of these scenarios, vector-specific transmission efficiency due to co-evolved begomovirus-vector complexes, has been implicated as a considerable driving force behind disease incidence, prevalence, and even severity (Idris et al. 2001, Brown and Czosnek 2002, Chowda-Reddy et al. 2012, Götz et al. 2012, Pan et al. 2012, Ning et al. 2015), and so, shifts in vector mitotype composition, particularly of exotics such as the B or Q mitotypes, can lead to changes in the prevalence of endemic viruses (Costa et al. 1993, Brown 2007). In addition, the invasive B and Q mitotypes have been found to be ‘promiscuous vectors’, transmitting begomoviruses at high efficiency that are endemic to the locales they have invaded, as well as those with which they have co-evolved (Davino et al. 2006, see refs in Brown 2007, Ghanim 2014, Orfanidou et al. 2016). Also, the host range preferences of different mitotypes can vary from monophagous to oligophagous, or polyphagous (Bird 1957, Mound 1963, Brown and Bird 1992, Bedford et al. 1994, Fishpool and Burban 1994, Brown et al. 1995b). A recent study of the Asia II-7 mitotype endemic to China Qiu et al. (2011) showed that development time, survivorship, sex ratio, and fecundity varied by plant host plant, i.e., collard, cucumber, hibiscus, laurel, poinsettia, and tomato (Qiu et al. 2011), underscoring the importance of understanding whitefly-plant-virus inter-relationships to enhance predictions of impending outbreaks, reduce the likelihood that outbreaks will occur, and better manage them when miscalculations have occurred. Comparison Between Fragments of the mtCOI Despite the extensive knowledge about the biology of B. tabaci, the systematics and the basis for delimitating species within the B. tabaci s.s.g. taxon have not been resolved. In most studies to date, the 3′-end of the COI gene (at ~780 bp) has been the molecular marker of choice for phylogenetic inference (Brown et al. 1995b, Brown 2010, Boykin et al. 2012, Liu et al. 2012, Hadjistylli et al. 2016), and by which unexpectedly high diversity has been predicted among B. tabaci mitotypes. However, no systematic comparison of the phylogenetic signal harbored by the 3′-fragment has been carried out until now. The tree topologies reconstructed herein using the 3′- and 5′-fragment sequences were overall in agreement with previous studies (Brown 2010, Dinsdale et al. 2010, Boykin et al. 2012). However, concatenating the 5′- and 3′-fragments to produce a marker of nearly twice the length, resulted in a more robust phylogenetic tree and with higher bv than either of the trees reconstructed for the partial COI gene fragments. Thus, among the three markers, the concatenated fragment spanning nearly the entire COI gene, at 1269 bases, produced the most robust phylogenetic tree. The results indicate that tree topologies for the 5′- and 3′-COI fragments for mitotypes from Pakistan were consistent at the level of major Asia I and Asia II clades. However, the sister clade relationships differed for the mitotypes that grouped in the Asia II and Sub-Saharan Africa major clades, depending on sequence diversity and number of sequences analyzed, with the divergence being most variable between the 3′- and 5′-ends for mitotypes grouped in the Asia II-7 sister clades. This observation is inconsistent with the conclusion from a previous study, which reported that the two fragments were equally informative (Ashfaq et al. 2014). However, it is in agreement with another that found the 3′-fragment could better resolved intraspecific relationships (Thomas and Ramamurthy 2014). And although the nt composition of the 5-′ and 3′-fragments was similar, as was also observed by Hsieh et al. (2014), the number of segregating sites, which is considered a useful measure of informativeness, was higher for the 3′-fragment. However, indicators of nt variability, including the haplotype diversity estimates based on the third codon positions, π diversity, and θ diversity at the first and second codon positions, were significantly higher for the 5′- than the 3′-fragment. Taken together the results of this study provide a convincing argument that the 3′-COI sequence, compared to the 5′-fragment, provides a more reliable gene tree for B. tabaci particularly at the mitotype level. The informativeness of the 3′ was further supported by evidence of extensive differences in the ratio of transitions over transversions, at 2.9 and 5.9, for the 3′ and 5′ fragment, respectively, which confirmed that transversions are the most abundant type of substitution in the 3′-COI fragment. It seems unlikely that the between-ratio differences can be accounted for by amino acid sequences dissimilarities, given that the ratio of degenerate sites at the 0-, two, and four-fold levels are similar among the two fragments (Table 5). This pattern is similar to the observations of Erpenbeck et al. (2006) who reported a Ti/Tv ratio for the 5′-COI fragment that was two times that of the 3′-fragment, and Lehr et al. (2005) found evidence for a high Ti/Tv ratio (6.8) in mtCOI sequences of the Anopheles albitarsus species complex. Conclusion The identification and genetic diversity studies of B. tabaci mitotypes can be accomplished using either the 3′- or 5′-COI fragment, or the concatenated sequence consisting of both. This study focused primarily on haplotypes that form three mitotype groups, extant in Pakistan affiliated with the Asia II major clade. Among them, the extent of variability in sequence polymorphisms and extent of divergence between the 3′- and 5′-COI sequences was reflected by weak phylogenetic support that was evident below the ‘major clade’ level, regardless of whether the COI fragment or the concatenated sequence was considered in the analysis. Because overall, the most robust phylogenietic support was obtained by analyzing an almost complete COI sequence (concatenated sequence), analysis of the entire COI gene sequence will better resolve relationships among and between both sister and major clades, and infer more reliable global COI gene trees, however, the trees cannot be reliably rooted. Second, this study provided robust evidence for recent genetic expansion of the Asia II-1 haplotypes, with widespread dispersion and establishment throughout the Punjab and Sindh provinces, establishing its predominance in the cropping systems. This phenomenon appears to have occurred at about the same time that CLCuKoV-Bur, a recombinant begomovirus, was documented to have emerged in the Punjab province, and spread throughout Pakistan and north India. These different lines of evidence strongly predict the recent demographic expansion of the polyphagous B. tabaci mitotypes, preceded by a genetic bottleneck. This bottleneck, resulting in an initial decrease in overall B. tabaci mitotype diversity, and reduced genetic variability in Asia II-1, appears to be directly associated with mitotype selection attributed to the overuse of pesticides to reduce the spread of highly damaging begomoviruses causing leaf curl disease in cotton and vegetable crops throughout Pakistan. The Asia II-7 mitotype on the other hand, appears to have been displaced in cotton-growing areas by the Asia II-1, while remaining a predominant type in the northern Punjab region near Lahore, where pesticide use is minimal and vegetation consists primarily of urban landscape and small-holding vegetable gardens where selection for insecticide resistance is probably rare or minimal. Acknowledgments The first author is supported on a fellowship provided by SENESCYT, (Ecuador) scholarship program 2011. This research was funded in part through the Pakistan-U.S. Cotton Productivity Enhancement Program of the USDA Agricultural Research Service, Project No.58-6402-0-178F. References Cited Ahmad , M. , and R. A. Khan . 2017 . Field-evolved resistance of Bemisia tabaci (Hemiptera: Aleyrodidae) to carbodiimide and neonicotinoids in Pakistan . J. Econ. Entomol . 110 : 1235 – 1242 . Google Scholar Crossref Search ADS PubMed Ahmad , M. , M. I. Arif , Z. Ahmad , and I. Denholm . 2002 . Cotton whitefly (Bemisia tabaci) resistance to organophosphate and pyrethroid insecticides in Pakistan . Pest Manag. Sci . 58 : 203 – 208 . Google Scholar Crossref Search ADS PubMed Ahmad , M. , M. I. Arif , and M. Naveed . 2010 . Dynamics of resistance to organophosphate and carbamate insecticides in the cotton whitefly Bemisia tabaci (Hemiptera: Aleyrodidae) from Pakistan . J. Pest Sci . 83 : 409 – 420 . Google Scholar Crossref Search ADS Ahmed , M. Z. , S. X. Ren , N. S. Mandour , M. N. Maruthi , M. Naveed , and B. L. Qiu . 2010 . Phylogenetic analysis of Bemisia tabaci (Hemiptera: Aleyrodidae) populations from cotton plants in Pakistan, China, and Egypt . J. Pest Sci . 83 : 135 – 141 . Google Scholar Crossref Search ADS Ahmed , M. Z. , P. J. De Barro , J. M. Greeff , S. X. Ren , M. Naveed , and B. L. Qiu . 2011 . Genetic identity of the Bemisia tabaci species complex and association with high cotton leaf curl disease (CLCuD) incidence in Pakistan . Pest Manag. Sci . 67 : 307 – 317 . Google Scholar Crossref Search ADS PubMed Alemandri , V. , P. De Barro , N. Bejerman , E. B. Argüello Caro , A. D. Dumón , M. F. Mattio , S. M. Rodriguez , and G. Truoli . 2012 . Species within the Bemisia tabaci (Hemiptera: Aleyrodidae) complex in soybean and bean crops in Argentina . J. Econ. Entomol . 105 : 48 – 53 . Google Scholar Crossref Search ADS PubMed Amrao , L. , S. Akhter , M. N. Tahir , I. Amin , R. W. Briddon , and S. Mansoor . 2010 . Cotton leaf curl disease in Sindh province of Pakistan is associated with recombinant begomovirus components . Virus Res . 153 : 161 – 165 . Google Scholar Crossref Search ADS PubMed Ashfaq , M. , P. D. Hebert , M. S. Mirza , A. M. Khan , S. Mansoor , G. S. Shah , and Y. Zafar . 2014 . DNA barcoding of Bemisia tabaci complex (Hemiptera: Aleyrodidae) reveals southerly expansion of the dominant whitefly species on cotton in Pakistan . PLoS One . 9 : e104485 . Google Scholar Crossref Search ADS PubMed Bandelt , H. J. , P. Forster , and A. Röhl . 1999 . Median-joining networks for inferring intraspecific phylogenies . Mol. Biol. Evol . 16 : 37 – 48 . Google Scholar Crossref Search ADS PubMed Banks , G. K. , J. Colvin , R. V. Chowda , and M. N. Maruthi . 2001 . First Report of the Bemisia tabaci B Biotype in India and an associated Tomato leaf curl virus disease epidemic . Plant Dis . 85 : 231 . Google Scholar Crossref Search ADS PubMed Bedford , I. D. , R. W. Briddon , J. K. Brown , R. C. Rosell , and P. G. Markham . 1994 . Geminivirus transmission and biological characterisation of Bemisia tabaci (Gennadius) biotypes from different geographic regions . Ann. Appl. Biol . 125 : 311 – 325 . Google Scholar Crossref Search ADS Benson , D. A. , M. Cavanaugh , K. Clark , I. Karsch-Mizrachi , D. J. Lipman , J. Ostell , and E. W. Sayers . 2012 . GenBank . Nucleic Acids Res . 41 : D36 – D42 . Google Scholar Crossref Search ADS PubMed Berry , S. D. , V. N. Fondong , C. Rey , D. Rogan , C. M. Fauquet , and J. K. Brown . 2004 . Molecular evidence for five distinct Bemisia tabaci (Homoptera: Aleyrodidae) geographic haplotypes associated with cassava plants in sub-Saharan Africa . Ann. Entomol. Soc. Am . 97 : 852 – 859 . Google Scholar Crossref Search ADS Bird , J. 1957 . A whitefly-transmitted mosaic of Jatropha gossypifolia . Agricultural Experimental Station. University of Puerto Rico. Tech. Paper . 22 : 1 – 35 . Boykin , L. M. , R. G. Shatters , Jr , R. C. Rosell , C. L. McKenzie , R. A. Bagnall , P. De Barro , and D. R. Frohlich . 2007 . Global relationships of Bemisia tabaci (Hemiptera: Aleyrodidae) revealed using Bayesian analysis of mitochondrial COI DNA sequences . Mol. Phylogenet. Evol . 44 : 1306 – 1319 . Google Scholar Crossref Search ADS PubMed Boykin , L. M. , K. F. Armstrong , L. Kubatko , and P. De Barro . 2012 . Species delimitation and global biosecurity . Evol. Bioinform. Online . 8 : 1 – 37 . Google Scholar Crossref Search ADS PubMed Briddon , R. W. , and P. G. Markham . 2000 . Cotton leaf curl virus disease . Virus Res . 71 : 151 – 159 . Google Scholar Crossref Search ADS PubMed Brower , A. 1994 . Phylogeny of Heliconius butterflies inferred from mitochondrial DNA sequences (Lepidoptera: Nymphalidae) . Mol. Phylogenet Evol . 3 : 159 – 174 . Google Scholar Crossref Search ADS PubMed Brown , J. K. 2000 . Molecular markers for the identification and global tracking of whitefly vector-Begomovirus complexes . Virus Res . 71 : 233 – 260 . Google Scholar Crossref Search ADS PubMed Brown , J. K. 2007 . The Bemisia tabaci complex: genetic and phenotypic variation and relevance to TYLCV-vector interactions , pp. 25 – 56 . In H. Czosnek (eds.), Tomato yellow leaf curl virus disease . Springer , Dordrecht, Netherlands . Brown , J. K. 2010 . Phylogenetic biology of the Bemisia tabaci sibling species group, pp. 31 – 67 . In P. A. Stansly and S. E. Naranjo (eds.), Bemisia: bionomics and management of a global pest . Springer , New York, NY . Brown , J. K. 2017 . National Plant Disease Recovery System Recovery Plan: Cotton leaf curl virus complex . USDA-ARS Office of Pest Management Policy National Plant Disease Recovery System . http://www.ars.usda.gov/research/docs.htm?docid=14271 (2013; revised 2017). Brown , J. K. , and J. Bird . 1992 . Whitefly-transmitted geminiviruses in the Americas and the Caribbean Basin: past and present . Plant Dis . 76 : 220 – 225 . Google Scholar Crossref Search ADS Brown , J. K. , and H. Czosnek . 2002 . Whitefly transmission of plant viruses . Adv. Bot. Res . 36 : 65 – 100 . Google Scholar Crossref Search ADS Brown , J. K. , S. A. Coats , I. D. Bedford , P. G. Markham , J. Bird , and D. R. Frohlich . 1995a . Characterization and distribution of esterase electromorphs in the whitefly, Bemisia tabaci (Genn.) (Homoptera: Aleyrodidae) . Biochem. Genet . 33 : 205 – 214 . Google Scholar Crossref Search ADS Brown , J. K. , D. R. Frohlich , and R. C. Rosell . 1995b . The sweetpotato or silverleaf whiteflies: biotypes of Bemisia tabaci or a species complex? Annu. Rev. Entomol . 40 : 511 – 534 . Google Scholar Crossref Search ADS Burban , C. , L. D. C. Fishpool , C. Fauquet , D. Fargette , and J. C. Thouvenel . 1992 . Host-associated biotypes within West African populations of the whitefly Bemisia tabaci (Genn), (Homoptera: Aleyrodidae) . J. Appl. Entomol . 113 : 416 – 423 . Google Scholar Crossref Search ADS Byrne , D. N. , and T. S. Bellows , Jr . 1991 . Whitefly biology . Annu. Rev. Entomol . 36 : 431 – 457 . Google Scholar Crossref Search ADS Cahill , M. , F. J. Byrne , K. Gorman , I. Denholm , and A. L. Devonshire . 1995 . Pyrethroid and organophosphate resistance in the tobacco whitefly Bemisia tabaci (Homoptera: Aleyrodidae) . Bull. Entomol. Res . 85 : 181 – 187 . Google Scholar Crossref Search ADS Chowda-Reddy , R. V. , M. Kirankumar , S. E. Seal , V. Muniyappa , and G. Valand . 2012 . Bemisia tabaci phylogenetic groups in India and the relative transmission efficacy of Tomato leaf curl Bangalore virus by an indigenous and an exotic population . J. Integr. Agric . 11 : 235 – 248 . Google Scholar Crossref Search ADS Chu , D. , F. H. Wan , Y. J. Zhang , and J. K. Brown . 2010 . Change in the biotype composition of Bemisia tabaci in Shandong Province of China from 2005 to 2008 . Environ. Entomol . 39 : 1028 – 1036 . Google Scholar Crossref Search ADS PubMed Cock , M. J. W. 1993 . Bemisia tabaci, an update 1986–1992 on the cotton whitefly with an annotated bibliography . CAB International Institute of Biol. Control , Ascot, United Kingdom. 78 pp. Costa , H. S. , and J. K. Brown . 1991 . Variation in biological characteristics and esterase patterns among populations of Bemisia tabaci and the association of one population with silverleaf symptom production . Entomol.Exp.Appl . 61 : 211 – 219 . Google Scholar Crossref Search ADS Costa , A. S. , and R. C. Russell . 1975 . Failure of Bemisia tabaci to breed on cassava plants in Brazil (Homoptera: Aleyrodidae) . Ciencia e Cultura Saõ Paulo . 2 : 388 – 390 . Costa , H. S. , J. K. Brown , and D. N. Byrne . 1991 . Life history traits of the whitefly, Bemisia tabaci (Homoptera: Aleyrodidae) on six virus-infected or healthy plant species . Environ. Entomol . 20 : 1102 – 1107 . Google Scholar Crossref Search ADS Costa , H. S. , Brown , J. K. , Sivasupramaniam , S. , and Bird , J. 1993 . Regional distribution, insecticide resistance, and reciprocal crosses between the ‘A’ and ‘B’ biotypes of Bemisia tabaci . Insect Sci. Applic . 14 : 127 – 138 . Darriba , D. , G. L. Taboada , R. Doallo , and D. Posada . 2012 . jModelTest 2: more models, new heuristics and parallel computing . Nat. Methods . 9 : 772 . Google Scholar Crossref Search ADS PubMed Davino , S. , C. Napoli , M. Davino , and G. P. Accotto . 2006 . Spread of Tomato yellow leaf curl virus in Sicily: partial displacement of another geminivirus originally present . Eur. J. Plant Pathol . 114 : 293 – 299 . Google Scholar Crossref Search ADS De Barro , P. J. , and P. J. Hart . 2000 . Mating interactions between two biotypes of the whitefly, Bemisia tabaci (Hemiptera: Aleyrodidae) in Australia . Bull. Entomol. Res . 90 : 103 – 112 . Google Scholar Crossref Search ADS PubMed De Barro , P. J. , J. W. Trueman , and D. R. Frohlich . 2005 . Bemisia argentifolii is a race of B. tabaci (Hemiptera: Aleyrodidae): the molecular genetic differentiation of B. tabaci populations around the world . Bull. Entomol. Res . 95 : 193 – 203 . Google Scholar Crossref Search ADS PubMed Delatte , H. , B. Reynaud , M. Granier , L. Thornary , J. M. Lett , R. Goldbach , and M. Peterschmitt . 2005 . A new silverleaf-inducing biotype Ms of Bemisia tabaci (Hemiptera: Aleyrodidae) indigenous to the islands of the south-west Indian Ocean . Bull. Entomol. Res . 95 : 29 – 35 . Google Scholar Crossref Search ADS PubMed Dinsdale , A. , L. Cook , C. Riginos , Y. M. Buckley , and P. De Barro . 2010 . Refined global analysis of Bemisia tabaci (Hemiptera: Sternorrhyncha: Aleyrododidae: Aleyrodidae) mitochondrial cytochrome oxidase 1 to identify species level genetic boundaries . Ann. Entomol. Soc. Am . 103 : 196 – 208 . Google Scholar Crossref Search ADS Dittrich , V. 1987 . Resistance and hormoligosis as driving forces behind pest outbreaks. pp. 169 – 181 . In K. Y. Brent and R. K. Alein (eds.), Rational pesticide use . Cambridge University Press , Cambridge, United Kingdom . Edgar , R. C. 2004 . MUSCLE: multiple sequence alignment with high accuracy and high throughput . Nucleic Acids Res . 32 : 1792 – 1797 . Google Scholar Crossref Search ADS PubMed Erpenbeck , D. , J. N. A. Hooper , and G. Wörheide . 2006 . CO1 phylogenies in diploblasts and the ‘Barcoding of life’ - are we sequencing a suboptimal partition? Mol. Ecol. Notes 6 : 550 – 553 . Google Scholar Crossref Search ADS Esterhuizen , L. L. , K. G. Mabasa , S. W. Van Heerden , H. Czosnek , J. K. Brown , H. Van Heerden , and M. E. Rey . 2013 . Genetic identification of members of the Bemisia tabaci cryptic species complex from South Africa reveals native and introduced haplotypes . J. Appl. Entomol . 137 : 122 – 135 . Google Scholar Crossref Search ADS Farooq , J. , A. Farooq , M. Riaz , M. R. Shahid , F. Saeed , M. S. Iqbal , T. Hussain , A. Batool , and A. Mahmood . 2014 . Cotton leaf curl virus disease a principle cause of decline in cotton productivity in Pakistan (a mini review) . Can. J. Plant Prot . 2 : 9 – 16 . Firdaus , S. , B. Vosman , N. Hidayati , E. D. Jaya Supena , R. G. Visser , and A. W. van Heusden . 2013 . The Bemisia tabaci species complex: additions from different parts of the world . Insect Sci . 20 : 723 – 733 . Google Scholar Crossref Search ADS PubMed Fishpool , L. D. C. , and C. Burban . 1994 . Bemisia tabaci: the whitefly vector of African cassava mosaic geminivirus . Trop. Sci . 34 : 55 – 72 . Folmer , O. , M. Black , W. Hoeh , R. Lutz , and R. Vrijenhoek . 1994 . DNA primers for amplification of mitochondrial cytochrome c oxidase subunit I from diverse metazoan invertebrates . Mol. March Biol. Biotechnol . 3 : 294 – 299 . Frohlich , D. R. , I. Torres-Jerez , I , I. D. Bedford , P. G. Markham , and J. K. Brown . 1999 . A phylogeographical analysis of the Bemisia tabaci species complex based on mitochondrial DNA markers . Mol. Ecol . 8 : 1683 – 1691 . Google Scholar Crossref Search ADS PubMed Fu , Y. X. 1997 . Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection . Genetics . 147 : 915 – 925 . Google Scholar PubMed Gawel , N. J. , and A. C. Bartlett . 1993 . Characterization of differences between whiteflies using RAPD-PCR . Insect Mol. Biol . 2 : 33 – 38 . Google Scholar Crossref Search ADS PubMed Ghanim , M. 2014 . A review of the mechanisms and components that determine the transmission efficiency of Tomato yellow leaf curl virus (Geminiviridae; Begomovirus) by its whitefly vector . Virus Res . 186 : 47 – 54 . Google Scholar Crossref Search ADS PubMed Gill , R. J. , and J. K. Brown . 2010 . Systematics of Bemisia and Bemisia relatives: can molecular techniques solve the Bemisia tabaci complex conundrum – A taxonomist’s viewpoint, pp. 5 – 29 . In P. A. Stansly and S. E. Naranjo (eds.), Bemisia: bionomics and management of a global pest . Springer , New York, NY . Götz , M. , S. Popovski , M. Kollenberg , R. Gorovitz , J. K. Brown , J. M. Cicero , H. Czosnek , S. Winter , and M. Ghanim . 2012 . Implication of Bemisia tabaci heat shock protein 70 in begomovirus-whitefly interactions . J. Virol . 86 ( 24 ): 13241 – 13252 . Google Scholar Crossref Search ADS PubMed Greathead , A. H. 1986 . Host plants , pp. 17 – 26 . In: Cock , M.J.W (eds.), Bemisia tabaci-A Literature survey on the cotton whitefly with an annotated bibliography . CAB International Institutes, Biological Control . Silwood Park, UK . Hadjistylli , M. , J. K. Brown , and G. K. Roderick . 2010 . Tools and recent progress in studying gene flow and population genetics of the Bemisia tabaci sibling species group, pp. 69 – 103 . In P. A. Stansly and S. E. Naranjo (eds.), Bemisia: bionomics and management of a global pest . Springer , New York, NY . Hadjistylli , M. , G. K. Roderick , and J. K. Brown . 2016 . Global population structure of a worldwide pest and virus vector: genetic diversity and population history of the Bemisia tabaci sibling species group . PLoS One . 11 : e0165105 . Google Scholar Crossref Search ADS PubMed Hajibabaei , M. , D. H. Janzen , J. M. Burns , W. Hallwachs , and P. D. Hebert . 2006 . DNA barcodes distinguish species of tropical Lepidoptera . Proc. Natl. Acad. Sci. U. S. A . 103 : 968 – 971 . Google Scholar Crossref Search ADS PubMed Hameed , S. , S. Hameed , M. Sadia , and S. A. Malik . 2012 . Genetic diversity analysis of Bemisia tabaci populations in Pakistan using RAPD markers . Electron. J. Biotech . 15 ( 6 ): 6 – 6 . Google Scholar Crossref Search ADS Harpending , H. C. 1994 . Signature of ancient population growth in a low-resolution mitochondrial DNA mismatch distribution . Hum. Biol . 66 : 591 – 600 . Google Scholar PubMed Hasegawa , M. , H. Kishino , and T. Yano . 1985 . Dating of the human-ape splitting by a molecular clock of mitochondrial DNA . J. Mol. Evol . 22 : 160 – 174 . Google Scholar Crossref Search ADS PubMed Hassan , I. , L. Amin , S. Mansoor , and R. W. Briddon . 2017 . Further changes in the cotton leaf curl disease complex: an Indication of Things to Come? Virus Genes . 53 ( 6 ): 759 – 761 . Google Scholar Crossref Search ADS PubMed Hebert , P. D. , M. Y. Stoeckle , T. S. Zemlak , and C. M. Francis . 2004 . Identification of Birds through DNA Barcodes . PLoS Biol . 2 : e312 . Google Scholar Crossref Search ADS PubMed Hsieh , C. H. , C. C. Ko , C. H. Chung , and H. Y. Wang . 2014 . Multilocus approach to clarify species status and the divergence history of the Bemisia tabaci (Hemiptera: Aleyrodidae) species complex . Mol. Phylogenet. Evol . 76 : 172 – 180 . Google Scholar Crossref Search ADS PubMed Hu , J. , P. De Barro , H. Zhao , J. Wang , F. Nardi , and S. S. Liu . 2011 . An extensive field survey combined with a phylogenetic analysis reveals rapid and widespread invasion of two alien whiteflies in China . PLoS One . 6 : e16061 . Google Scholar Crossref Search ADS PubMed Huelsenbeck , J. P. , and F. Ronquist . 2001 . MRBAYES: Bayesian inference of phylogenetic trees . Bioinformatics . 17 : 754 – 755 . Google Scholar Crossref Search ADS PubMed Hurst , G. D. , and F. M. Jiggins . 2005 . Problems with mitochondrial DNA as a marker in population, phylogeographic and phylogenetic studies: the effects of inherited symbionts . Proc. Biol. Sci . 272 : 1525 – 1534 . Google Scholar Crossref Search ADS PubMed Ibáñez , C. M. , Cubillos , L. A. , Tafur , R. , Argüelles , J. , Yamashiro , C. and E. Poulin . 2011 . Genetic diversity and demographic history of Dosidicus gigas (Cephalopoda: Ommastrephidae) in the Humboldt Current System . March Ecol. Prog. Ser . 431 : 163 – 171 . Google Scholar Crossref Search ADS Idris , A. M. , Smith , S. E. , and Brown , J. K. 2001 . Ingestion, transmission, and persistence of Chino del tomate virus (CdTV), a New World begomovirus, by Old and New World biotypes of the whitefly vector Bemisia tabaci . Ann. Appl. Biol . 139 : 145 – 154 . Google Scholar Crossref Search ADS Islam , W. , W. Lin , M. Qasim , S. U. Islam , H. Ali , M. Adnan , M. Arif , Z. Du , and Z. Wu . 2018 . A nation-wide genetic survey revealed a complex population structure of Bemisia tabaci in Pakistan . Acta Trop . 183 : 119 – 125 . Google Scholar Crossref Search ADS PubMed Ivanova , N. V. , T. S. Zemlak , R. H. Hanner , and P. D. Hebert . 2007 . Universal primer cocktails for fish DNA barcoding . Mol. Ecol. Notes 7 : 544 – 548 . Google Scholar Crossref Search ADS Jones , D. R. 2003 . Plant viruses transmission by whiteflies . Eur. J. Plant Pathol . 109 : 195 – 219 . Google Scholar Crossref Search ADS Kimura , M. 1980 . A simple method for estimating evolutionary rates of base substitutions through comparative studies of nucleotide sequences . J. Mol. Evol . 16 : 111 – 120 . Google Scholar Crossref Search ADS PubMed Kumar , S. , G. Stecher , and K. Tamura . 2016 . MEGA7: molecular evolutionary genetics analysis version 7.0 for bigger datasets . Mol. Biol. Evol . 33 : 1870 – 1874 . Google Scholar Crossref Search ADS PubMed Lee , W. , J. Park , G. S. Lee , S. Lee , and S. Akimoto . 2013 . Taxonomic status of the Bemisia tabaci complex (Hemiptera: Aleyrodidae) and reassessment of the number of its constituent species . PLoS One . 8 : e63817 . Google Scholar Crossref Search ADS PubMed Legg , J. P. , P. Sseruwagi , S. Boniface , G. Okao-Okuja , R. Shirima , S. Bigirimana , G. Gashaka , H. W. Herrmann , S. Jeremiah , H. Obiero , et al. 2014 . Spatio-temporal patterns of genetic change amongst populations of cassava Bemisia tabaci whiteflies driving virus pandemics in East and Central Africa . Virus Res . 186 : 61 – 75 . Google Scholar Crossref Search ADS PubMed Legg , J. P. , R. French , D. Rogan , G. Okao-Okuja , and J. K. Brown . 2002 . A distinct Bemisia tabaci (Gennadius) (Hemiptera: Sternorrhyncha: Aleyrodidae) genotype cluster is associated with the epidemic of severe cassava mosaic virus disease in Uganda . Mol. Ecol . 11 : 1219 – 1229 . Google Scholar Crossref Search ADS PubMed Lehr , M. A. , C. W. Kilpatrick , R. C. Wilkerson , and J. E. Conn . 2005 . Cryptic Species in the Anopheles (Nyssorhynchus) albitarsis (Diptera: Culicidae) Complex: incongruence Between Random Amplified Polymorphic DNA-Polymerase Chain Reaction Identification and Analysis of Mitochondrial DNA COI Gene Sequences . Ann. Entomol. Soc. Am . 98 : 908 – 917 . Google Scholar Crossref Search ADS PubMed Li , S. J. , X. Xue , M. Z. Ahmed , S. X. Ren , Y. Z. Du , J. H. Wu , A. G. Cuthbertson , and B. L. Qiu . 2011 . Host plants and natural enemies of Bemisia tabaci (Hemiptera: Aleyrodidae) in China . Insect Sci . 18 ( 1 ): 101 – 120 . Google Scholar Crossref Search ADS Librado , P. , and J. Rozas . 2009 . DnaSP v5: a software for comprehensive analysis of DNA polymorphism data . Bioinformatics . 25 : 1451 – 1452 . Google Scholar Crossref Search ADS PubMed Liu , S. S. , P. J. De Barro , J. Xu , J. B. Luan , L. S. Zang , Y. M. Ruan , and F. H. Wan . 2007 . Asymmetric mating interactions drive widespread invasion and displacement in a whitefly . Science . 318 : 1769 – 1772 . Google Scholar Crossref Search ADS PubMed Liu , S. , J. Colvin , and P. J. De Barro . 2012 . Species concepts as applied to the whitefly Bemisia tabaci systematics: how many species are there? J. Integr. Agric . 11 ( 2 ): 176 – 186 . Google Scholar Crossref Search ADS Lopez-Avila , A. 1986 . Taxonomy and biology, pp. 3 . In M. J. W. Cock (ed.), Bemisia tabaci - A literature survey on the cotton whitefly with an annotated bibliography . CAB International Institute of Biological Control , Ascot, United Kingdom . Maddison , W. P. , and D. R. Maddison . 2014 . Mesquite: a modular system for evolutionary analysis, version 2.75. http://mesquiteproject.org Mahmood , T. , M. Arshad , M. I. Gill , H. T. Mahmood , M. Tahir , and S. Hussain . 2003 . Burewala strain of cotton leaf curl virus: a threat to CLCuV cotton resistant varieties . Asian J. Plant Sci . 2 : 968 – 970 . Google Scholar Crossref Search ADS Mansoor , S. , I. Amin , S. Iram , M. Hussain , Y. Zafar , K. A. Malik , and R. W. Briddon . 2003a . Breakdown of resistance in cotton to cotton leaf curl disease in Pakistan . Plant Pathol . 52 : 784 Google Scholar Crossref Search ADS Mansoor , S. , R. W. Briddon , S. E. Bull , I. D. Bedford , A. Bashir , M. Hussain , M. Saeed , M. Y. Zafar , K. A. Malik , C. Fauquet , and P. G. Markham . 2003b . Cotton leaf curl disease is associated with multiple monopartite begomoviruses supported by single DNA β . Arch. Virol . 148 : 1969 – 1986 . Google Scholar Crossref Search ADS Martin , J. H. , and L. A. Mound . 2007 . An annotated checklist of the world’s whiteflies (Insecta: Hemiptera: Aleyrodidae) . Zootaxa . 149 : 21 – 84 . Masood , M. , I. Amin , I. Hassan , S. Mansoor , J. K. Brown , and R. W. Briddon . 2017 . Diversity and distribution of cryptic species of the Bemisia tabaci (Hemiptera: Aleyrodidae) complex in Pakistan . J. Econ. Entomol . 20 : 1 – 6 . Masters , B. C. , V. Fan , and H. A. Ross . 2011 . Species Delimitation–a Geneious plugin for the exploration of species boundaries . Mol. Ecol. Resour . 11 : 154 – 157 . Google Scholar Crossref Search ADS PubMed Mirol , P. M. , J. Routtu , A. Hoikkala , and R. K. Butlin . 2008 . Signals of demographic expansion in Drosophila virilis . BMC Evol. Biol . 8 : 1 . Google Scholar Crossref Search ADS PubMed Mohanty , A. K. , and A. N. Basu . 1986 . Effects of host plants and seasonal factors on intraspecific variations in pupal morphology of the whitefly vector Bemisia tabaci (Genn.) (Homoptera: Aleyrodidae) . J. Entomol. Res . 10 : 19 – 26 . Mound , L. A. 1963 . Host correlated variation in Bemisia tabaci (Gennadius) (Homoptera:Aleyrodidae) . Proc. R. Entomol. Soc. London . 38 : 171 – 180 . Mound , L. A. , and S. H. Halsey . 1978 . Whitefly of the World: a systematic catalogue of the Aleyrodidae (Homoptera) with host plant and natural enemy data . British Museum and John Wiley & Sons , Chichester, United Kingdom . 340 pp. Moya , A. , P. Guirao , D. Cifuentes , F. Beitia , and J. L. Cenis . 2001 . Genetic diversity of Iberian populations of Bemisia tabaci (Hemiptera: Aleyrodidae) based on random amplified polymorphic DNA-polymerase chain reaction . Mol. Ecol . 10 : 891 – 897 . Google Scholar Crossref Search ADS PubMed Naveen , N. C. , R. Chaubey , D. Kumar , K. B. Rebijith , R. Rajagopal , B. Subrahmanyam , and S. Subramanian . 2017 . Insecticide resistance status in the whitefly, Bemisia tabaci genetic groups Asia-I, Asia-II-1 and Asia-II-7 on the Indian subcontinent . Sci. Rep . 7 : 40634 . Google Scholar Crossref Search ADS PubMed Ning , W. , X. Shi , B. Liu , H. Pan , W. Wei , Y. Zeng , X. Sun , W. Xie , S. Wang , Q. Wu , et al. 2015 . Transmission of Tomato Yellow Leaf Curl Virus by Bemisia tabaci as Affected by Whitefly Sex and Biotype . Sci. Rep . 5 : 10744 . Google Scholar Crossref Search ADS PubMed Orfanidou , C. G. , P. G. Pappi , K. E. Efthimiou , N. I. Katis , and V. I. Maliogka . 2016 . Transmission of Tomato chlorosis virus (ToCV) by Bemisia tabaci biotype Q and evaluation of four weed species as viral sources . Plant Dis . 100 : 2043 – 2049 . Google Scholar Crossref Search ADS PubMed Pan , H. , D. Chu , W. Yan , Q. Su , B. Liu , S. Wang , Q. Wu , W. Xie , X. Jiao , R. Li , et al. 2012 . Rapid spread of tomato yellow leaf curl virus in China is aided differentially by two invasive whiteflies . PLoS One . 7 : e34817 . Google Scholar Crossref Search ADS PubMed Paradis , E. , J. Claude , and K. Strimmer . 2004 . APE: analyses of Phylogenetics and Evolution in R language . Bioinformatics . 20 : 289 – 290 . Google Scholar Crossref Search ADS PubMed Perring , T. M. , A. D. Cooper , R. J. Rodriguez , C. A. Farrar , and T. S. Bellows , Jr . 1993 . Identification of a whitefly species by genomic and behavioral studies . Science . 259 : 74 – 77 . Google Scholar Crossref Search ADS PubMed Pons , J. , T. G. Barraclough , J. Gomez-Zurita , A. Cardoso , D. P. Duran , S. Hazell , S. Kamoun , W. D. Sumlin , and A. P. Vogler . 2006 . Sequence-based species delimitation for the DNA taxonomy of undescribed insects . Syst. Biol . 55 : 595 – 609 . Google Scholar Crossref Search ADS PubMed Qiu , B. L. , S. A. Coats , S. X. Ren , A. M. Idris , X. U. Caixia , and J. K. Brown . 2007 . Phylogenetic relationships of native and introduced Bemisia tabaci (Homoptera: Aleyrodidae) from China and India based on mtCOI sequencing and host plant comparisons . Prog. Nat. Sci . 17 : 645 – 654 . Google Scholar Crossref Search ADS Qiu , B. L. , F. Dang , S. J. Li , M. Z. Ahmed , F. L. Jin , S. X. Ren , and A. G. Cuthbertson . 2011 . Comparison of biological parameters between the invasive B biotype and a new defined Cv biotype of Bemisia tabaci (Hemiptera: Aleyradidae) in China . J. Pest. Sci . 84 : 419 – 427 . Google Scholar Crossref Search ADS Rajagopalan , P. A. , A. Naik , P. Katturi , M. Kurulekar , R. S. Kankanallu , and R. Anandalakshmi . 2012 . Dominance of resistance-breaking cotton leaf curl Burewala virus (CLCuBuV) in northwestern India . Arch. Virol . 157 : 855 – 868 . Google Scholar Crossref Search ADS PubMed Rambaut , A. , M. A. Suchard , D. Xie , and A. J. Drummond . 2014 . Tracer v1.6. http://beast.bio.ed.ac.uk/Tracer Ramírez-Soriano , A. , S. E. Ramos-Onsins , J. Rozas , F. Calafell , and A. Navarro . 2008 . Statistical power analysis of neutrality tests under demographic expansions, contractions and bottlenecks with recombination . Genetics . 179 : 555 – 567 . Google Scholar Crossref Search ADS PubMed Ramos-Onsins , S. E. , and J. Rozas . 2002 . Statistical properties of new neutrality tests against population growth . Mol. Biol. Evol . 19 : 2092 – 2100 . Google Scholar Crossref Search ADS PubMed Rogers , A. R. , and H. Harpending . 1992 . Population growth makes waves in the distribution of pairwise genetic differences . Mol. Biol. Evol . 9 : 552 – 569 . Google Scholar PubMed Ronquist , F. , M. Teslenko , P. van der Mark , D. L. Ayres , A. Darling , S. Höhna , B. Larget , L. Liu , M. A. Suchard , and J. P. Huelsenbeck . 2012 . MrBayes 3.2: efficient Bayesian phylogenetic inference and model choice across a large model space . Syst. Biol . 61 : 539 – 542 . Google Scholar Crossref Search ADS PubMed Rosell , R. C. , I. D. Bedford , D. R. Frohlich , R. J. Gill , J. K. Brown , and P. G. Markham . 1997 . Analysis of morphological variation in distinct populations of Bemisia tabaci (Homoptera: Aleyrodidae) . Ann. Entomol. Soc. Am . 90 : 575 – 589 . Google Scholar Crossref Search ADS Ross , H. A. , S. Murugan , and W. L. Li . 2008 . Testing the reliability of genetic methods of species identification via simulation . Syst. Biol . 57 : 216 – 230 . Google Scholar Crossref Search ADS PubMed Rozas , J. , A. Ferrer-Mata , J. C. Sánchez-DelBarrio , S. Guirao-Rico , P. Librado , S. E. Ramos-Onsins , and A. Sánchez-Gracia . 2017 . DnaSP 6: DNA sequence polymorphism analysis of large data sets . Mol. Biol. Evol . 34 : 3299 – 3302 . Google Scholar Crossref Search ADS PubMed Ruegg , K. C. , and T. B. Smith . 2002 . Not as the crow flies: a historical explanation for circuitous migration in Swainson’s thrush (Catharus ustulatus) . Proc. Biol. Sci . 269 : 1375 – 1381 . Google Scholar Crossref Search ADS PubMed Russell , L. M. 1948 . The North American species of whiteflies of the genus Trialeurodes . Misc. Publ. U.S. Dep. Agric. No. 635, Washington DC , 85 pp. Russell , L. M. 1957 . Synonyms of Bemisia tabaci (Gennadius) (Homoptera, Aleyrodidae) . Bull. Brooklyn Entomol. Soc . 52 : 122 – 123 . Sattar , M. N. , Z. Iqbal , M. N. Tahir , and S. Ullah . 2017 . The Prediction of a New CLCuD Epidemic in the Old World . Front. Microbiol . 8 : 631 . Google Scholar Crossref Search ADS PubMed Shimodaira , H. , and M. Hasegawa . 1999 . Multiple comparisons of log-likelihoods with applications to phylogenetic inference . Mol. Biol. Evol . 16 : 1114 – 1116 . Google Scholar Crossref Search ADS Simon , C. , F. Frati , A. Beckenbach , B. Crespi , H. Liu , and P. Flook . 1994 . Evolution, weighting, and phylogenetic utility of mitochondrial gene sequences and a compilation of conserved polymerase chain reaction primers . Ann. Entomol. Soc. Am . 87 : 651 – 701 . Google Scholar Crossref Search ADS Simón , B. , J. L. Cenis , F. Beitia , S. Khalid , I. M. Moreno , A. Fraile , and F. García-Arenal . 2003 . Genetic structure of field populations of begomoviruses and of their vector Bemisia tabaci in Pakistan . Phytopathology . 93 : 1422 – 1429 . Google Scholar Crossref Search ADS PubMed Singer , M. S. 2008 . Evolutionary ecology of polyphagy, pp. 29 – 42 . In K. J. Tilmon (ed.), Specialization, speciation and radiation: evolutionary biology of herbivorous insects . University of California Press , Berkeley, CA . Slatkin , M. , and R. R. Hudson . 1991 . Pairwise comparisons of mitochondrial DNA sequences in stable and exponentially growing populations . Genetics . 129 : 555 – 562 . Google Scholar PubMed Smith , M. A. , B. L. Fisher , and P. D. Hebert . 2005 . DNA barcoding for effective biodiversity assessment of a hyperdiverse arthropod group: the ants of Madagascar . Philos. Trans. R. Soc. Lond. B. Biol. Sci . 360 : 1825 – 1834 . Google Scholar Crossref Search ADS PubMed Song , H. , J. E. Buhay , M. F. Whiting , and K. A. Crandall . 2008 . Many species in one: DNA barcoding overestimates the number of species when nuclear mitochondrial pseudogenes are coamplified . Proc. Natl. Acad. Sci. U. S. A . 105 : 13486 – 13491 . Google Scholar Crossref Search ADS PubMed Sseruwagi , P. , J. Legg , M. N. Maruthi , J. Colvin , M. E. C. Rey , and J. K. Brown . 2005 . Genetic diversity of Bemisia tabaci (Gennadius) (Hemiptera: Aleyrodidae) populations and presence of the B biotype and a non-B biotype that can induce silverleaf symptoms in squash, in Uganda . Ann. Appl. Biol . 147 : 253 – 265 . Google Scholar Crossref Search ADS Sseruwagi , P. , M. N. Maruthi , J. Colvin , M. E. C. Rey , J. K. Brown , and J. P. Legg . 2006 . Colonization of non-cassava plant species by cassava whiteflies (Bemisia tabaci) in Uganda . Entomol. Exp. Appl . 119 : 145 – 153 . Google Scholar Crossref Search ADS Swofford , D. L. 2003 . PAUP*. Phylogenetic Analysis Using Parsimony (*and Other Methods). Version 4 . Sinauer Associates , Sunderland, MA . Tajima , F. 1989 . Statistical method for testing the neutral mutation hypothesis by DNA polymorphism . Genetics . 123 : 585 – 595 . Google Scholar PubMed Takahashi , R. 1936 . Some Aleyrodidae, Aphididae, Coccidae (Homoptera), and Thysanoptera from Micronesia . Tenthredo . 1 : 109 – 120 . Tariq , M. I. , S. Afzal , I. Hussain , and N. Sultana . 2007 . Pesticides exposure in Pakistan: a review . Environ. Int . 33 : 1107 – 1122 . Google Scholar Crossref Search ADS PubMed Tay , W. T. , S. Elfekih , L. N. Court , K. H. J. Gordon , H. Delatte , and P. J. De Barro . 2017 . The trouble with MEAM2: implications of pseudogenes on species delimitation in the globally invasive Bemisia tabaci (Hemiptera: Aleyrodidae) cryptic species complex . Genome Biol. Evol . 9 : 2732 – 2738 . Google Scholar Crossref Search ADS PubMed Thomas , A. , and V. V. Ramamurthy . 2014 . Multiple gene markers to understand genetic diversity in the Bemisia tabaci (Hemiptera: Aleyrodidae) species complex . Florida Entomol . 97 : 1451 – 1457 Google Scholar Crossref Search ADS Villesen , P. 2007 . FaBox: an online toolbox for fasta sequences . Mol. Ecol. Notes 7 : 965 – 968 . Google Scholar Crossref Search ADS Viscarret , M. M. , I. Torres-Jerez , E. Agostini de Manero , S. N. López , E. E. Botto , and J. K. Brown . 2003 . Mitochondrial DNA evidence for a distinct clade of New World Bemisia tabaci (Genn.) (Hemiptera: Aleyrodidae) from Argentina and Bolivia, and presence of the Old World B biotype in Argentina . Ann. Entomol. Soc. Am . 96 : 65 – 72 . Google Scholar Crossref Search ADS Ward , R. D. , T. S. Zemlak , B. H. Innes , P. R. Last , and P. D. Hebert . 2005 . DNA barcoding Australia’s fish species . Philos. Trans. R. Soc. Lond. B. Biol. Sci . 360 : 1847 – 1857 . Google Scholar Crossref Search ADS PubMed Wickham , H. 2009 . ggplot2: Elegant Graphics for Data Analysis . Springer-Verlag , New York, NY . Zaidi , S. S. , M. Shafiq , I. Amin , B. E. Scheffler , J. A. Scheffler , R. W. Briddon , and S. Mansoor . 2016 . Frequent occurrence of Tomato Leaf Curl New Delhi Virus in cotton leaf curl disease affected cotton in Pakistan . PLoS One . 11 : e0155520 . Google Scholar Crossref Search ADS PubMed Zhang , Y.-P. , J. K. Uyemoto , and BC Kirkpatrick . 1998 . A small scale procedure for extracting nucleic acids from woody plants infected with various phytopathogens for PCR assay . J. Virol. Methods 71 : 45 – 50 . Google Scholar Crossref Search ADS PubMed Zwickl , D. J. 2006 . Genetic algorithm approaches for the phylogenetic analysis of large biological sequence datasets under the maximum likelihood criterion. Ph.D. dissertation . The University of Texas at Austin , Austin, TX . © The Author(s) 2019. Published by Oxford University Press on behalf of Entomological Society of America. All rights reserved. For permissions, please e-mail: journals.permissions@oup.com. This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/journals/pages/open_access/funder_policies/chorus/standard_publication_model)
TI - Demographic Expansion of the Predominant Bemisia tabaci (Gennadius) (Hemiptera: Aleyrodidae) Mitotypes Associated With the Cotton Leaf Curl Virus Epidemic in Pakistan
JF - Annals of the Entomological Society of America
DO - 10.1093/aesa/saz002
DA - 2019-05-07
UR - https://www.deepdyve.com/lp/oxford-university-press/demographic-expansion-of-the-predominant-bemisia-tabaci-gennadius-6KogP0oRwR
SP - 265
VL - 112
IS - 3
DP - DeepDyve
ER -