The Bemisia tabaci (Gennadius) (Hemiptera, Aleyrodidae) sibling (or cryptic) species group is comprised of an unknown number of morphologically indistinguishable phenotypic variants, for which a small number have been characterized biologically (Bedford et al. 1994; Brown et al. 1995; Sseruwagi et al. 2005, 2006; Caballero 2007; Gill & Brown 2010; Hadjistylli et al. 2016). As a species group, they are endemic to tropical, subtropical, and temperate ecosystems. As a group, they also have an extensive host range among eudicots (Mound & Halsey 1978; Cock 1993), with more than 200 plant species recorded as hosts in Pakistan (Attique et al. 2003).

Recently, B. tabaci mitotypes and their genetic groups have been differentiated based on phylogenetic analysis of an approximately 580 to 780 base pair fragment of the 3′-end of the mitochondrial cytochrome oxidase I gene (mtCOI-??3′) sequence (Frohlich et al. 1999), with several theoretical estimates proposed for species delineation at 3.5 to 4.0% nucleotide (nt) divergence, irrespective of corresponding biological species criteria (Dinsdale et al. 2010; Lee et al. 2013). To date, no consensus has been reached for a globally applicable species demarcation system within the B. tabaci taxon. Even so, understanding the mitotype composition and population structure of B. tabaci in agroecosystems is important, because management practices may cause shifts in the prevalence of the different B. tabaci mitotypes, which in turn may govern begomovirus transmission competency, potentially influencing the dynamics of the different mitotypes and cotton leaf curl disease outbreaks, possibly resulting in shifting of 1 or more components of the pathosystem. Different mitotypes of the B. tabaci sibling species group have been shown to respond differently to various pesticide regimes, the range of plants suitable as hosts, and environmental factors. These and other factors strongly influence the composition of B. tabaci mitotypes and their dynamics in agroecosystems (Bedford et al. 1994; Brown 2007, 2010; Gill & Brown 2010).

Plant viruses transmitted by whiteflies belong to 1 of 5 genera (Jones 2003; Navas-Castillo et al. 2011). Among the economically important plant viruses transmitted by whiteflies, the genus Begomovirus (Geminiviridae) is the most widespread, comprising more than 400 species ( https://talk.ictvonline.org/taxonomy/), including many known to cause economically important plant diseases (Brown et al. 2012). In Pakistan and India, over 60 begomoviral species or strains have been described from cultivated and wild plant species (Ho et al. 2017). Cotton leaf curl disease is caused by one of several monopartite begomoviruses or strains, and associated betasatellites that are endemic to much of southern Asia, including Pakistan (Zhou et al. 1998; Saleem et al. 2016). During the past several decades, the emergence of cotton leaf curl disease has resulted in major economic losses in the cotton producing areas of Pakistan (Briddon et al. 2001; Zubair et al. 2017). Understanding the distribution of the different whitefly mitotypes and their population dynamics, which can be predicted by spatio-temporal patterns of mtCOI diversity, have become essential for devising effective management approaches. Such strategies are further enhanced by knowledge of the host range, pesticide resistance, and other biological characteristics associated with the range of different mitotypes that transmit the equally diverse begomoviral species and strains that cause cotton leaf curl disease, as well as mitotype-virus specific factors that directly or indirectly influence the efficiency of whitefly-mediated transmission.

The objective of this study was to determine the mitotype composition and spatial patterns of genetic divergence of the B. tabaci whitefly vector collected from commercial cotton fields in the southern Punjab districts of Pakistan from 2014 to 2016. Mitotypes of the whitefly B. tabaci were identified phylogenetically based on the mtCOI-3′ gene sequence, and the structure of haplotypes was studied using a median-joining network analysis.

Materials and Methods

WHITEFLY DNA ISOLATION, POLYMERASE CHAIN REACTION AMPLIFICATION, AND DNA SEQUENCING

Total DNA was purified from individual whiteflies using a previously published protocol (Zhang et al. 1998) with modifications (Paredes-Montero et al. 2019). Briefly, whiteflies were dipped in double-distilled water to dilute the ethanol, followed by blotting against the torn edge of a piece of No. 1 Whatman filter paper (Sigma-Aldrich, Darmstadt, Germany) to remove the distilled water. Each whitefly was ground with a micro-pestle in the bottom of a 1.0 mL microfuge tube, containing 250 µL lysis buffer. The lysis buffer contained 100 mM Tris-HCl, 1.4 M NaCl, 20 mM EDTA, 2% hexadecyltrimethylammonium bromide (CTAB) (Sigma-Aldrich, Darmstadt, Germany), and 20 µg per reaction of proteinase K. The lysate was incubated at 55 °C for 1 h followed by 10 min at 65 °C. One volume of chloroform:isoamyl-alcohol (24:1) was added, the contents were mixed by inversion, and the emulsion was separated by centrifugation, at 14,000 rpm for 3 min. The DNA was precipitated by the addition of 1 volume isopropanol and 40 µg glycogen, followed by incubation at 4 °C for 10 min. The pellet was collected by micro-centrifugation at 14,000 rpm for 10 min, washed with 70% ethanol, air dried, dissolved in 30 µL ddH₂O, and stored at –20 °C.

To detect the possible presence of the exotic B mitotype [North Africa-Mediterranean-Middle East major clade [NAF-MED-ME (III)] among the endemic Asia I and II mitotypes (major clades) (Brown 2010), samples were subjected to polymerase chain reaction amplification using the previously reported B-mitotype-specific primers 5′-BtBF1F-3′ [TATTTCACTTCAGCCACTATAA] and 5′-WfBr2R-3′ [GCTTAAATCTTACTAACCGCAG], yielding an expected size fragment of 550 base pairs (Andreason et al. 2017). The negative samples were subjected to polymerase chain reaction amplification using the universal mtCOI primers C1-J-2195F 5′-TTGATTTTTTGGTCATCCAGAAGT-3′ and L2-N-3014R 5′-TCCAATGCACTAATCTGCCATATTA-3′ that yield an expected size amplicon of about 850 base pairs for B. tabaci (Frohlich et al. 1999). The polymerase chain reaction cycling parameters were 95 °C for 2 min, followed by 30 cycles at 95 °C for 60 s, 52 °C for 60 s, 72 °C for 60 s, with a final extension at 72 °C for 5 min. Reactions contained 1X Jumpstart REDTaq Ready-Mix (Sigma-Aldrich, St. Louis, Missouri, USA), 0.4 µM primers, 20 ng DNA, and ddH₂O, to a final volume of 25 µL.

Fig. 1.

Map of the Punjab Province in Pakistan showing the locations of whitefly sample collections, and distribution of the Asia II-1 mitotype of Bemisia tabaci, represented by the star symbol.

The amplicons were ligated into the pGEM®-T Easy plasmid vector (Promega Corp., Madison, Wisconsin, USA). The presence of an expected size insert of 875 to 880 base pairs (non-B mitotype) and 550 base pairs (B mitotype), respectively, was confirmed by colony polymerase chain reaction amplification (Gussow & Clackson 1989). Amplicons were sequenced bi-directionally by capillary Sanger DNA sequencing (University of Arizona Genetics Core, University of Arizona, Tucson, Arizona, USA). The sequences were assembled, trimmed, manually edited, and aligned using Lasergene software (DNASTAR, Madison, Wisconsin, USA).

Results

NETWORK ANALYSIS FOR ASIA II-1 HAPLOTYPES

Based on the median-joining network, 4 predominant haplotypes were resolved within the Punjab region Asia II-1 population, and the haplotypes in cotton crops were found to be ‘geographically unrestricted' or independent of district. Also, several rare haplotypes were identified, based on the accumulation of unique mutations, relative to the ‘common’ haplotype. The rare haplotypes in the network were identified based on their affiliation with a star-like pattern surrounding the common haplotypes (Fig. 4). This kind of pattern can be indicative of population recovery, which is characteristically reflected as increased genetic diversity following a decline in diversity within a population (Fig. 4).

Fig. 2.

Phylogenetic tree showing the relationship among members of the Bemisia tabaci complex found in Asia, based on the 3′-mitochondria cytochrome oxidase I gene (mtCOI-3′). All whiteflies collected for this study during 2014 to 2016, were identified as the Asia II-1 mitotype. The sequence of the Bemisia afer species was used as ‘outgroup’ to root the tree. The sequences in the upper-most collapsed clade clustered all of the sequences obtained in this study.

Discussion

The results indicated that Asia II-1 group sub-mitotype 1 was the only B. tabaci variant among the samples collected from the cotton-growing areas of the Punjab Province. This observation is consistent with recent studies that reported it to be the predominant mitotype in agroecosystems in the Punjab and Sindh provinces (Ashfaq et al. 2014; Masood et al. 2017; Islam et al. 2018; Paredes-Montero et al. 2019). Five other endemic mitotypes or sub-mitotypes have been documented previously in Pakistan agroecosystems, namely Asia I and Asia II-1, Asia II-5, Asia II-7, and Asia II-8 (Firdaus et al. 2013; Masood et al. 2017; Islam et al. 2018; Paredes-Montero et al. 2019). Additionally, there is one report of the B mitotype (Ashfaq et al. 2014) which is phylogeographically affiliated with the NAF-MED-ME (III) major clade, and so is thought to have been introduced only recently (Ashfaq et al. 2014).

The diversity of haplotypes given by the corrected pairwise distances varied by district. The most highly divergent haplotypes at 1.37% were found in Khanewal, whereas in Dera Ghazi Khan, haplotypes differed only by 0.5%. The range of genetic divergence among whiteflies from the different districts studied was found to be similar to divergences reported in several previous studies (Ashfaq et al. 2014; Masood et al. 2017; Paredes-Montero et al. 2019), suggesting that differences in whitefly management practices may have favored or limited establishment or displacement of certain haplotypes compared to others (Chu et al. 2010; Horowitz & Ishaaya 2014). Certain classes of insecticides over others have been shown to favor outbreaks of different B. tabaci mitotypes/sub-mitotypes, which appears to be due to differences in their ability to develop insecticide resistance (Bedford et al. 1994; Coats et al. 1994; Anthony et al. 1995; Denholm et al. 1998; Ahmad et al. 2002; Dennehy et al. 2010; Luo et al. 2010; Horowitz & Ishaaya 2014).

Fig. 3.

Corrected pairwise distances and number of haplotypes within Asia II-1 sequences per district. The error bars indicate the extent of the standard error.

The network of haplotypes exhibited a high frequency of unique mutations that gave rise to several low-frequency haplotypes that diverged from the predominant haplotype in 1 or 2 mutations, forming a star-like structure (Mirol et al. 2008). This pattern is suggestive of genetic expansion or recovery from a recent bottleneck within the Asia II-1 populations. A previous study conducted in these same locations from 2012 to 2014 showed a similar scenario, with haplotype diversity levels starting to recover after a dramatic decline in genetic diversity among Asia II-1 haplotypes starting in the early 2000s (Paredes-Montero et al. 2019). It was hypothesized that this pattern of reduced genetic variation could have been due to the reasonably well-documented overuse of pesticides in cotton agroecosystems over the last 2 decades (Ahmad et al. 2002). Two yr after initial signs of recovery occurred, the rate of genetic expansion was about the same for B. tabaci collected from 2014 to 2016; however, the number of predominant haplotypes increased from 1 to 4 in relation to that reported by Paredes-Montero et al. (2019). In this context, and despite the smaller sample size analyzed in this study, at 107 compared to 538 samples collected in many of the same or additional districts analyzed by Paredes-Montero et al. (2019), the genetic signatures indicated the continued predominance of haplotypes classified in the subclade Asia II-1 in cotton grown in the Punjab Province from 2014 to 2016. Therefore, based on these observations, sustained genetic expansion of the populations is predicted to reach a plateau in subsequent yr in the absence of changes in cotton varieties or management practices.

Fig. 4.

Median-joining network analysis of Asia II-1 sequences (= 0). Yellow circles indicate Bemisia tabaci haplotypes, and the size of the circle reflects relative abundance. The black lines indicate the extents of relatedness between haplotypes, and the red numbers placed on each line indicate the number of mutations.

The exclusive presence of Asia II-1 in cotton-agroecosystems in the southern Punjab Province may explain in part the recent spread of several predominant strains and species of cotton leaf curl disease there (Paredes-Montero et al. 2019). This hypothesis is supported by the high incidence of cotton leaf curl disease observed in cotton plantings where the Asia II-1 mitotype predominates. In contrast, lower disease incidence has been reported in the southern province of Sindh (Simón et al. 2003; Ahmed et al. 2011; Pan et al. 2018), where the B biotype occurs (Ashfaq et al. 2014) together with mitotypes endemic to the Asian continent (Ashfaq et al. 2014; Masood et al. 2017; Islam et al. 2018; Paredes-Montero et al. 2019). Based on one study, the (exotic) B mitotype showed low efficiency to no transmission of an isolate of the cotton leaf curl multan virus introduced recently in China (Pan et al. 2018), compared to variable rates of transmission by the introduced Q mitotype, or high rates of transmission by mitotypes endemic to Asia. This evidence for differential transmission specificity, together with a low transmission frequency of cotton leaf curl multan virus-transmission by the B mitotype suggests that, although the introduced B mitoype has established in at least 1 locale in Pakistan, it is unlikely to have a prominent role in the transmission of the majority of begomoviruses occurring in Pakistan. However, the B miotype has been shown to be a competent vector of other New and Old World begomoviruses, and so it could pose a threat as a begomovirus vector if any of those viruses are introduced into Pakistan in the future.

Conclusion

The Asia II-1 was identified as the predominant mitotype among B. tabaci samples collected from cotton in southern Punjab, Pakistan. Here, 2 lines of evidence underscore differences in diversification among and between different B. tabaci mitotypes associated with the different districts. First, the intra-mitotype diversity was found to vary by district, which most likely is indicative of differential mitotype responsiveness to locally imposed agricultural practices. Second, a continued increase in haplotype diversity, or genetic expansion, was observed following a notable decline in diversity beginning in about 2010 with the initiation of recovery during 2014 to 2016 that could in time lead to restored diversity, perhaps approaching the previously high levels, and the geographic redistribution of the Asia II-1 and other sub-mitotypes in the southern Punjab Province of Pakistan. Several interacting factors are thought to be involved, including altered management practices implemented since 2010, including use of different pesticide regimes for whitefly control, particularly those that cause minimal agroecosystem and natural enemy disturbances. Additionally, with the release of tolerant or resistant cotton leaf curl disease cotton varieties, transmission frequencies may have been reduced, which if so encourages the routinely judicious use of insecticides. Routine monitoring of the population structure and distribution of B. tabaci mitotypes and cotton leaf curl disease-begomoviruses in cotton or cotton-vegetable agroecosystems in the Punjab Province is expected to contribute to implementing best practices that contribute to avoiding upsurges in whitefly that lead to cotton leaf curl disease outbreaks. Consequently, the new information presented here as well as that from previous studies of historical fluxes in B. tabaci mitotype and begomovirus distribution are extremely valuable for predicting these outbreaks.