Analysis of the bacterial communities and endosymbionts of natural populations of Bemisia tabaci in several crop fields from Mexico semi-arid zone

Abstract

Purpose

Bemisia tabaci (Aleyrodidae family) is an insect vector of plant viruses that affects a wide variety of crops around the world. In the following study, we analyzed the variation in the bacterial communities associated with natural populations of B. tabaci (MEAM1) of four different crops from six regions from Mexico semi-arid zone (Baja California Sur).

Methods

PCR was used to amplify the mitochondrial cytochrome oxidase I gene (mtCOI), and then to carry out the phylogenetic analysis for genetic identification of the isolated B. tabaci. Next generation sequencing coupled with 16S metagenomic analysis was applied in order to characterize B. tabaci inner microbial community. Finally, bacterial obligate symbiont and facultative symbiont were confirmed by PCR amplification and by phylogenetic analysis.

Results

Ours results pointed toward that B. tabaci MEAM1 inner bacterial communities were predominantly structured by Proteobacteria phylum. Moreover, the most represented endosymbionts were the obligate endosymbiont from the genus “Candidatus Portiera”, as well as two facultative symbionts belonging to genera Rickettsia and Hamiltonella; both obligate and facultative endosymbionts were present for all samples, and their relative abundance varied was crop-independent.

Conclusions

Geographic localization and insect diet play a central role to maintain bacterial community structure of the B. tabaci MEAM1 whitefly at phylum taxonomic level. Agricultural practices were a factor that affected samples of bacterial community structure similarities, reflected in samples clustering. Host plants that are part of B. tabaci diet did not influence directly, in spite of sap nutrient differences, into obligate and facultative endosymbionts relative abundance.

Introduction

Animals are a nutrient-rich habitat for microorganisms with the ability to modulate animal immune system and metabolic adaptations to colonize and utilize host resources (Kostic et al. 2013; Chung et al. 2017). The association between insects and microbial organisms is quite important for ecological diversification and evolution of insects’ feeding habits. The next generation sequencing (NGS) has emerged as an important tool to describe the microbiota structure from environmental and biological samples such as insect pest, suggesting that associated microbiota plays an important role in protection against plants defenses, resistance to insecticides, temperature, reproduction, and even contributing in host nutrition (Weiss and Aksoy 2011; Su et al. 2016; Mereghetti et al. 2017; Chung et al. 2017).

The ability of insects to adapt into certain nutritional lifestyles can be attributed to survival strategies that also include associations with obligate symbionts or “primary bacterial symbionts” which have been characterized in several insect groups, such as Buchnera in aphids, “Candidatus Portiera” in whiteflies and Carsonella in psyllids. The associations of insects with facultative endosymbionts or “secondary endosymbionts” are not essential for host survival, such as Rickettsia, Hamiltonella, and Cardinium in whiteflies and Arsenophonus in psyllids and aphids, but these endosymbionts can manipulate hosts in many other ways, such as genetic differentiation, adaptation to a wide range of food plants, defense against natural enemies, reproduction, insecticide resistance, and pathogen transmission (Chuche et al. 2017).

Bemisia tabaci is recognized as a cryptic species complex, a small sap-sucking insect known as whitefly that is capable of damaging a wide range of agricultural products through its feeding activity and its active role in plant viruses transmission (De Barro et al. 2011). Phylogenetic analyses applied for B. tabaci were based on the mitochondrial cytochrome oxidase I gene marker (mtCOI) and have been shown that there are four major geographically defined clades: (i) sub-Saharan African, (ii) New World, (iii) Asia, and (iv) Africa/Middle East/Asia Minor/Central Asia/Mediterranean. From those regions, the B. tabaci varieties that belong to clade IV, Mediterranean (MED), and Middle East-Asia Minor 1 (MEAM1) have become into some of the most damaging agricultural pests, which are biologically differentiated and reproductively isolated from each other (Boykin et al. 2007, 2013). Both B. tabaci biotype species (MED and MEAM1) have also been reported in the American continent, with frequent invasion into several crop fields giving rise to genetic material exchange process and increasing pest traits (Elfekih et al. 2018). B. tabaci is an insect that its mouthparts were developed for piercing and sucking, and derived by its feeding habits has established a symbiotic relationship with bacteria of genus “Candidatus Portiera”, that is an obligate symbiont, but also with facultative symbionts of several genera including Rickettsia, Hamiltonella, Wolbachia, Arsenophonus, Cardinium, and Fritschea (Gnankiné et al. 2013). Recent progress has been made in whitefly-plants interactions to determine the role of facultative symbionts, such as Hamiltonella that plays a central role abolishing tomato plant defense, and also has been related in Tomato yellow leaf curl virus (TYLCV) transmission (Su et al. 2013, 2015). In addition, it has been determined that Rickettsia presence is directly related to heat tolerance and protection against the entomopathogenic Pseudomonas syringae in B. tabaci (Brumin et al. 2011; Hendry et al. 2014).

Most of insect symbionts are vertically transmitted although recent works have provided evidence of horizontal transmission by experiments with axenic larvae that were fed with plants possessing Rickettsia, and then those insects in adult stages were analyzed by PCR to detect Rickettsia 16S rDNA gene sequence (Chrostek et al. 2017; Li et al. 2017a; Li et al. 2017b). The prevalence of facultative symbionts is often attributed to multiple factors, including climate, geographical localization, geographical isolation (surrounded by deserts, mountains, rivers, sea, absence of appropriate host plants), horizontal transmission, and insect host plant, in natural populations of B. tabaci, and other insects (Pan et al. 2012; Jing et al. 2014; Zhu et al. 2018). The insect microbial communities have been characterized mainly by PCR-based detection techniques focused on 16S rDNA or 23S rDNA gene amplification. Nowadays, the 16S rDNA metagenomic analysis of phloem sap from which insects were fed has opened the landscape about exploring ecological questions. Specifically, in B. tabaci, it has been determined that the associated bacterial community structure does not differ significantly between insects bred in laboratory and collected from crop fields, and, in both insect populations, obligate symbionts and facultative symbionts were the most abundant (Jing et al. 2014). Furthermore, the work with B. tabaci fed in laboratory with several plant species shed light into the essential role of host plant for shaping the whitefly-associated bacterial community structure (Su et al. 2016).

Baja California peninsula in Mexico offers climatological, geographical, and agricultural features that favor high populations of whitefly, such as relative geographical isolation and intensive agriculture in well-defined and isolated areas with a low diversity of plant species that favor the bacteria able to colonize the whitefly gut, and a semi-arid weather to study whiteflies’ natural populations and their associated bacterial communities. We analyzed the variation in the bacterial communities and endosymbionts associated with natural populations of B. tabaci of four different crops, from six regions, from this semi-arid zone (Baja California Sur) of Mexico to assess the impact of host plants, geographical localization, and agricultural practices (Fig. S1).

Materials and methods

Whitefly population collection

The whiteflies used in this study were natural populations of B. tabaci collected in six sites from the Baja California peninsula (south region), Mexico, during the time period of 2016–2017. The crop plant species (host-plant) localization were (i) La Matanza (23° 39′ 26.2″ N, 110° 16′ 00.3″ W) with Solanum lycopersicum (1T) and Capsicum annuum var. serrano (2CH) (two samples); (ii) El Comitán (24° 08′ 06.6″ N, 110° 25′ 40.3″ W) with C. chinense (4H) and C. annuum var. poblano (5P) (two samples); (iii) Los Planes (23° 55′ 58.1″ N, 109° 54′ 33.1″ W) with Cucurbita pepo (3CO) and Brassica oleracea var. gemmifera (8CA) (two samples); (iv) El Valle (25° 04′ 04.7″ N, 111° 36′ 17.4″ W) with C. pepo (two samples: 6CA and 7CA); (v) El Carrizal (23° 48′ 01.1″ N, 110° 17′ 39.1″ W) with C. pepo (one sample: 10CA); and (vi) El Pescadero (23° 23′ 39.8″ N, 110° 10′ 28.9″ W) with C. pepo (one sample: 9CA); totaling ten samples (Table 1). The whiteflies (in adult stage) collected were carried out during their feeding periods and placed in 95% (v/v) ethanol until total DNA extraction.

Table 1 Descriptive information of each sample collected in Baja California peninsula in Mexico used in this study

DNA extraction

Whiteflies in adult stage (white color insect body, red eyes) were previously identified with the aid of a stereoscope, for every localization and crop (Azab et al. 1969; Rosell et al. 1995). Twenty-one whiteflies were pooled per localization and host-plant (sample) (21 insects × 10 pool samples = 21 insects/pool sample = 210 insects). Those insects positively identified were processed and described as follows: each insect pooled sample was washed three times with 70% (v/v) ethanol and the supernatant was removed; each pooled sample was processed with plastic micropestle and 800 μL of DNA extraction buffer (2% CTAB, 1.4 M NaCl, 20 mM EDTA, 100 mM Tris-HCl pH 8.0, 0.2% β-mercaptoethanol, 1 mg/mL Proteinase K) (Chen et al. 2015), incubated at 55 °C for 60 min with occasional stirring, followed by the addition of one volume of chloroform:isoamyl alcohol (24:1), and centrifuged for 15 min at 12,000×g. After that, one volume of isopropyl alcohol (Sigma-Aldrich, USA) was added and incubated overnight at − 20 °C. DNA was precipitated by centrifugation at 12,000×g for 10 min, and then washed with 70% (v/v) ethanol followed by centrifugation at 12,000×g. The DNA pellet was resuspended in 40 μL of ultra-pure distilled water. Total DNA samples were treated with RNase A (10 mg/mL) (Promega, Madison, WI, USA) and incubated at 37 °C for 30 min. Total DNA (100–200 ng) integrity for each sample was verified by agarose 0.9% (w/v) gel electrophoresis, and its quality was assessed with the ratio of λ 260 nm/280 nm (from 1.79 to 2.00) with a spectrophotometer NanoDrop 2000 (Thermo Fisher Scientific, Waltham, MA, USA). Finally, the samples were stored at − 20 °C for its preservation.

Detection of B. tabaci species (mtCOI gene)

The mitochondrial cytochrome oxidase I (mtCOI) gene from whitefly samples was amplified using biotype-specific primers reported for Shatters et al. (2009). PCR reactions were performed using 0.25 mM of each primer, 0.5-mM deoxynucleotide triphosphates (dNTPs), 1.5-mM MgCl₂, 0.5 μL of Taq DNA polymerase buffer, 1 u of GoTaq (Promega, Madison, WI, USA), and 100 ng of DNA. PCR reaction was carried out as follows: initial denaturation step (95 °C, 5 min, 1 cycle), amplification step for 35 cycles (95 °C, 30 s; 64 °C, 1 min and 72 °C, 1 min), and a final extension at 72 °C for 5 min. PCR amplification was purified using the Zymoclean™ Gel DNA Recovery Kit (Zymo Research, Irvine, California, USA) and Sanger sequenced. The accession number for the B. tabaci identified is MK169413.

16S V3 rDNA sequencing

DNA of 21 whiteflies were pooled and processed as one sample that represents a localization and host-plant. A total of ten DNA pools were used to amplify the V3-V4 region of the 16S rDNA gene using the primers V3-338F and V4-806R. The product was used as template for the second round of amplification using the primers V3-338F and V3-533r containing Illumina adapters. Samples were sequenced on the Illumina MiniSeq platform 2 × 150 (300 cycles). The resulting sequences were deposited in the GenBank under the BioProject accession number PRJNA521383.

Metagenomic analysis

For relative abundance analysis of endosymbionts, unassembled DNA sequences were annotated with Metagenomics Rapid Annotation using Subsystems Technology (MG-RAST) pipeline version 3.5 for downstream analyses (htt://metagenomic.anl.gov; Meyer et al. 2008). Taxonomic assignments were performed using the normalized abundance of sequence matches to M5nr and SILVA SSU database. For data normalization, the frequency of representative hits to each individual taxon for each metagenome was divided by the total number of hits. To identify hits, Blastn was employed with a maximum E-value of 1e-05, a minimum identity of 97%, and maximum alignment of 15 bp obtained from MG-RAST server was downloaded and analyzed. The sequences were deposited in MG-RAST under the project “whitefly Baja California” (ID 82439).

In addition, Illumina MiniSeq sequencing reads were subjected to quality filtering (Q ≥ 33) using the modified Mott algorithm; sequences were paired-end, and then merged with BBmerge (Ewing et al. 1998; Ewing and Green 1998; Kearse et al. 2012; Bushnell et al. 2017). Those sequence reads were then de novo assembled with Geneious Assembler algorithm of Geneious R10.2 (www.geneious.com), with a minimum identity of 98%, and both contigs and unassembled reads were compared against the GenBank non-redundant (“nt-nr”) database using Megablast (Randle-Boggis et al. 2016). Megablast results were used to create a curated database specific for each pool sample. Sequence reads per pool sample were clustered against its corresponding created database using Sequence Classifier plugin, with a minimum of 99% identity for species taxonomic assignation, to finally obtain the table of operational taxonomic units (OTUs) frequencies using Geneious R10.2, (www.geneious.com). The resulting OTUs tables were then processed to be analyzed with R programming language using a variety of packages and custom scripts (www.r-project.org). We estimated Chao1, Shannon, and Simpson (alpha) biodiversity indices with the package “vegan”. PERMANOVA statistical analysis was performed with “adonis” function with the package “vegan”, and ANOVA analysis followed with Tukey post hoc test with “HSD.test” function with the package “agricolae”. Bray-Curtis distance estimations were calculated using the “vegdist” function, as well as principal component analysis using “prcomp” function with the package “ggfortify” (Oksanen et al. 2014; Coleman-Derr et al. 2016; Castañeda and Barbosa 2017).

Detection of whitefly symbionts

The obligate symbiont (bacteria from genus “Candidatus Portiera”) and facultative symbionts (bacteria from genera Rickettsia, Hamiltonella, Cardinium, Wolbachia, Fritschea, Arsenophonus) PCR-based detections were performed using specific primer sequences (Portiera: 28F and 1098R, Rickettsia: Rb-F and Rb-R, Hamiltonella: Ham-F and Ham-R, Cardinium: CFB-F and CFB-R, Wolbachia: wsp81-F and wsp691-R, Fritschea: U23F and 23SIGR, Arsenophonus: Ars23S-1 and Ars23S-2) reported in Pan et al. (2012). PCR reactions were performed with 100 ng of DNA as follows: initial denaturation step (95 °C, 5 min, 1 cycle), amplification step for 35 cycles (95 °C, 30 s denaturation; annealing 55 °C (Wolbachia)), 58 °C (“Candidatus Portiera”, Cardinium, Arsenophonus, and Hamiltonella), 60 °C (Rickettsia and Fritschea for 1 min; and extension 72 °C, 1 min), and a final extension at 72 °C for 5 min. The PCR amplification products were purified using the Zymoclean™ Gel DNA Recovery Kit (Irvine, California, USA) and Sanger sequenced. The accession numbers for 16S rDNA sequences of the assessed samples are Rickettsia (MK159392), Portiera (MK208462), “Candidatus Hamiltonella defensa” (MK280777).

Phylogenetic analysis of species of B. tabaci and bacterial endosymbionts

The mtCOI gene for B. tabaci, and 16S rDNA and 23S rDNA for bacterial endosymbionts sequences were obtained from GenBank and with the consensus sequence from each sample (three obtained sequences per sample were aligned and assembled to yield a consensus sequence) obtained in this work to be phylogenetically analyzed. The sequences alignment was carried out with MUSCLE with default parameters, and phylogenetic analysis was performed with maximum likelihood (ML) methodology coupled with Kimura-2 parameter model, and 1500 replicates of bootstrap for statistical support of evolutionary branch lengths, with Molecular Evolutionary Genetic Analysis software package version 7 (MEGA7) (Kumar et al. 2016).

Results

Genetic identification of B. tabaci based on mtCOI

All B. tabaci collected samples were addressed through PCR-based biotype species identification using the mtCOI gene marker, as well as their corresponding phylogenetic analysis to support the identification. The resulting PCR-based species identification for all B. tabaci samples were positive for biotype MEAM1 (Fig. S2). B. tabaci mtCOI marker gene phylogenetic analysis confirmed that whiteflies’ natural populations of the Baja California peninsula south region belong to biotype MEAM1 (MEAM1_MEX_BCS) by the monophyletic grouping with other biotype MEAM1, and also had an evolutionary close relationship with the B. tabaci MEAM1 isolate from Brazil (Fig. 1).

Fig. 1

Phylogenetic tree based on sequence of cytochrome oxidase I (mtCOI) fragment of B. tabaci collected in Baja California peninsula, Mexico (MEAM1_MEX_BCS), and selected reference sequences from GenBank. The alignment and phylogenetic analysis by MUSCLE alignments coupled with maximum-likelihood (ML) methodology, based on Kimura-2 parameter model. The length of branches is supported by 1500 bootstrap replicates. Trialeurodes vaporariorum was used as outgroup

Microbial community diversity and composition of B. tabaci MEAM1

To characterize the bacterial community structures in B. tabaci biotype MEAM1_MEX_BCS natural populations from Baja California peninsula south region (Table 1), we proceed with 16S rDNA metagenomic analysis. The analysis was conducted for ten samples which represent four different crops from six regions. The number of yielded raw reads per pool sample was as follows: 1T, 160686 reads; 2CH, 162,672 reads; 3CO, 179,140 reads; 4H, 165,424 reads; 5P, 180,306 reads; 6CA, 170,192 reads; 7CA, 147,540 reads; 8CA, 200,446 reads; 9CA, 194,206 reads; and 10CA, 119,364 reads. We obtained 755,360 paired-end sequence reads after quality filtering (from 1,679,976 total reads, 1,510,720 reads were above a phred of 33, and then paired-end merged), for all (ten) collected samples. The percentage of reads did not match to any known operational taxonomic unit (OTU) in the GenBank nucleotide non-redundant (“nt-nr”) database was of 2.602%. We were able to obtain 2093 OTUs, and after singleton removal, we obtained 1209 OTUs.

The alpha diversity indices were estimated to determine bacterial community abundance (observed OTUs, Chao1 index), richness (Shannon), and homogeneity (Simpson) for B. tabaci MEAM1 collected samples. Bacterial community abundance indicated by observed OTUs and Chao1 index estimation was relatively low, as expected, by the low OTUs abundance that were represented by a single read, for all B. tabaci MEAM1 collected samples, and supported by the Shannon index (2.533 ± 0.207). Interestingly, the Simpson index (0.94 ± 0.02) showed a homogenous bacterial composition for all B. tabaci MEAM1 collected samples (Table 2). It is worth noting that both, Chao1 index and Simpson index, reflect the stability of gut microbiota as a discrete and well-defined ecosystem to promote the proliferation of certain and specific bacterial groups, and only influenced by the insects’ feeding habits.

Table 2 Summary of sequencing analysis of 16S rDNA gene of B. tabaci MEAM1 and alpha diversity indices

Furthermore, we found for all collected samples ten phyla such as including Proteobacteria, Acidobacteria, Actinobacteria, Firmicutes, Bacteroidetes, Chloroflexi, Cyanobacteria, Planctomycetes, Tenericutes, and Verrucomicrobia. The dominant phylum in B. tabaci MEAM1_MEX_BCS was Proteobacteria (the relative abundances of 96–99%). Two classes of Proteobacteria phylum were the most abundant: Gammaproteobacteria (58–75%) and Alphaproteobacteria (0.060–40%). The Gammaproteobacteria class was mainly conformed by two of genera “Candidatus Portiera” (24–54%) and Hamiltonella (17.93–44%). The Alphaproteobacteria class was mainly conformed by the genus Rickettsia (0.011–39.34%,). The remaining phyla also exhibited low abundance frequencies (~ 3–< 0.6%) (Fig. 2). The endosymbionts (obligate and facultative) of B. tabaci MEAM1_MEX_BCS belong to Proteobacteria phylum, and this fact could represent a major contribution into crop production, since it is possible to develop strategies focused against those (endosymbionts) bacteria to decrease or eliminate their abundance in crop fields.

Fig. 2

Stacked column bar graph representing bacterial community composition of natural population of B. tabaci MEAM1_MEX_BCS based on the relative abundance at a phylum level and b family level (enlarged part of the low abundance associated bacteria)

We procced with β-diversity analysis based in Bray-Curtis distances estimation to determine sample hierarchical clustering, and a principal component analysis (PCA) to address the impact of host plant species, localizations, and agrochemicals application (agricultural practices) on bacterial community structure of B. tabaci MEAM1. PCA analysis showed three discrete and well-defined groups or samples clustering: (i) 10CA, 1T, 9CA, and 8CA; (ii) 7CA, 3CO, 4H, 6CA, and 5P; and (iii) 2CH (Fig. 3). To test the effect of host plant species, localization of crop fields, and agricultural practices on the observed variations on bacterial community structure of B. tabaci MEAM1, permutational multivariate analysis of variance (PERMANOVA) was performed. PERMANOVA confirmed that agrochemicals application was the main factor (host plant species and localization were not statistically significant) that exerted the effect on the bacterial community structure clustering (R² = 0.2971, P = 0.007), determined in PCA analysis (Fig. 3). These sample grouping was also supported by the hierarchical clustering analysis derived from the Bray-Curtis distances estimation (Fig. S3), and also by an analysis of correlation among samples and heatmap analysis (Fig. S4 and Fig. S5, respectively).

Fig. 3

Principal component analysis (PCA) of B. tabaci MEAM1_MEX_BCS bacterial community structures from natural populations of Baja California peninsula (south region), Mexico

Host plant effect on whitefly symbiont abundances

The whitefly symbionts relative abundance was assessed from eight samples collected, excepting both 3CO (family Brassicaceae) and 7CA (“El Valle” repeated sample) samples, and divided into two groups based on the plant host family: Solanaceae (1T, 2CH, 4H, 5P) and Cucurbitaceae (6CA, 8CA, 9CA, 10CA). The whitefly symbionts were identified at genus taxonomic level for “Candidatus Portiera”, Rickettsia, and Hamiltonella. To test the effect of host plant family (Solanaceae and Cucurbitaceae) on whitefly symbiont abundances, a two-way ANOVA with a Tukey post hoc test was performed. The most abundant whitefly symbiont was obligate symbiont “Candidatus Portiera” with a relative abundance range of 31.61 to 54.37%; showing a slight, but not significant, higher abundance in insects collected from Solanaceae (46.4–54.37%) than those collected from Cucurbitaceae (31.61–51.05%) (Fig. 4). The facultative symbionts showed similar abundance for both plant host families (Rickettsia, 0.011–33.69% from Solanaceae, and 22.46–34.80% from Cucurbitaceae; and Hamiltonella, 17.93–44.2% from Solanaceae, and 17.72–31.45% from Cucurbitaceae). The obligate symbiont (“Candidatus Portiera”) overall relative abundance was significantly higher than the facultative symbionts (Rickettsia and Hamiltonella) relative abundances (P = 0.0009); however, host plant family did not have an effect on both obligate and facultative symbionts abundances (P > 0.05) (Table 3 and Fig. 4).

Fig. 4

Comparison of relative abundance of the endosymbionts “Candidatus Portiera”, Hamiltonella, and Rickettsia between the B. tabaci MEAM1_MEX_BCS collected from Solanaceae and Cucurbitaceae

Table 3 Relative abundance (%) of bacterial symbionts in B. tabaci MEAM1

Detection and species identification of symbionts from B. tabaci MEAM1 natural populations

The bacterial symbionts (obligate and facultative) in B. tabaci MEAM1_MEX_BCS samples were also confirmed with PCR-based detection. We were able to detect the presence of the bacterial symbionts (obligate symbiont “Candidatus Portiera”, and facultative symbionts Rickettsia and Hamiltonella) for all samples analyzed. The PCR products were sequenced, and those sequences were analyzed for species identification with “blastn” program. The obligate symbiont “Candidatus Portiera aleyrodidarum” was identified (CP016304.1, 91% identity), as well as facultative symbionts Rickettsia sp. MEAM1 (CP016305.1), and “Candidatus Hamiltonella defense” (CP017613), with 98% and 94% identity, respectively. Phylogenetic analysis for the sequences corresponding for the three symbiont bacteria detected and identified in B. tabaci MEAM1_MEX_BCS for “Candidatus Portiera aleyrodidarum”, “Candidatus Hamiltonella defense”, and Rickettsia sp. MEAM1 supported the metagenomic analysis results (Figs. 5, 6, and 7).

Fig. 5

Phylogenetic tree based on sequence of 16S rDNA of “Candidatus Portiera” fragment obtained from B. tabaci MEAM1_MEX_BCS, and selected reference sequences from GenBank. The alignment and phylogenetic analysis by MUSCLE alignments coupled with maximum-likelihood (ML) methodology, based on Kimura-2 parameter model. The length of branches is supported by 1000 bootstrap replicates. Buchnera aphidicola endosymbiont of Schizaphis graminum was used as outgroup

Fig. 6

Phylogenetic tree based on sequence of 16S rDNA of Rickettsia endosymbiont fragment obtained from B. tabaci collected in Baja California peninsula, Mexico (MEAM1_MEX_BCS), and selected reference sequences from GenBank. The number on the branches represents bootstrap values for 1500 replicates. Wolbachia endosymbiont of B. tabaci was used as outgroup sequences

Fig. 7

Phylogenetic tree based on sequence of 16S rDNA of “Candidatus Hamiltonella” fragment obtained from B. tabaci MEAM1_MEX_BCS, and selected reference sequences from GenBank. The length of branches is supported by 1500 bootstrap replicates. “Candidatus Regiella” and “Candidatus Serratia” endosymbiont of aphids was used as outgroup

Discussion

The whitefly B. tabaci is an insect pest species for a wide range of crops comprising several biotypes, which are based in mtCOI gene marker. The B. tabaci feeding habit is almost restricted to plant sap obtained directly from phloem, and as a consequence, the insect requires a close association with its symbiotic bacteria. Several studies have been focused to understand the relationship of B. tabaci with its symbiotic bacteria, and also have been determined that there are multiple factors affecting this interaction (Pan et al. 2012; Su et al. 2016). In addition, a study conducted in a controlled environment, the diet influence was addressed onto B. tabaci bacterial community structures, as well as the insect biotype (Su et al. 2016). And also, it has been determined that insect geographic localization is another main factor, that impacts in the insect diet, to influence on the presence of the associated symbiotic bacteria (Pan et al. 2012). We took into account that the three factors mentioned above (insect biotype, diet, and geography) influence directly onto the B. tabaci bacterial community structure. Hence, in the present study, we characterized natural populations of whiteflies (B. tabaci) in the Baja California peninsula south region, Mexico, applying several approaches such as identification of insect biotype, analysis of bacterial community structures associated to insects, and the identification of the obligate symbiont and facultative symbiont from insect samples (Fig. S1).

We were able to collect whitefly samples of natural populations from several localizations in the Baja California peninsula south region, and we determined that those collected whiteflies belong to biotype MEAM1 (Fig. 1). It is worth noting that this result was consistent with previous studies which have been shown that since the introduction of B. tabaci biotype MEAM1 in the American continent in the late 1980s, the native B. tabaci biotype NW has been displaced. The B. tabaci biotype MEAM1 dominance is mainly attributed to traits such as insecticide resistance, high-level polyphagia, and fecundity (Schuster et al. 2010; Smith et al. 2016; MycKenzie and Osborne 2017).

Once we have determined that B. tabaci natural populations from Baja California peninsula south region belong to biotype MEAM1, we proceeded to characterize the bacterial community structure associated with the whitefly by 16S rDNA metagenomic analysis. The most abundant phylum was Proteobacteria (96–99%) in the bacterial community of B. tabaci MEAM1_MEX_BCS, the remaining phyla (nine) represented abundances from ~ 3% to less than 0.6% (Fig. 2), indicating a low richness (Shannon index) and a high homogeneity (Simpson index) (Table 2). Those results were quite similar with previous reports which were determined that phloem sap from different plant sources diet of insects contributes directly onto low bacterial richness and a homogenous structure (Jing et al. 2014; Zhao et al. 2016).

The bacterial community structures for all B. tabaci MEAM1 samples were mainly composed by the phylum Proteobacteria; then we proceed to analyze global differences and similarities among all samples, to determine whether some factors may influence on to their possible clustering or not. Interestingly, the B. tabaci MEAM1 samples were clustered in two well-defined groups, and just one sample apparently did not cluster. The sample clustering was influenced by agriculture practice for crop fields, since samples 1T, 8CA, 9CA, and 10 CA that came from crop fields with agrochemical (pesticides/organic pesticides) regular application were clustered (Fig. 3); and samples, 4H, 5P, 6CA, and 7CA that came from crop fields without agrochemical application along with 3CO of the Brassicaceae family were clustered in another group (Fig. 3). However, the sample 2CH, which corresponds to C. annuum from crop field with agrochemical application, was grouped apart from sample 4H and 5P, which also corresponds to C. annuum from crop field without agrochemical application, and also grouped apart from other two groups; probably, the agrochemical application induces a significant metabolic change in this plant species that also affected the B. tabaci MEAM1 bacterial community structure, as reflected in the PCA analysis (Fig. 3). This result was consistent with Zhao et al.’s (2016) study, which determined that the bacterial community structures associated to natural population of aphids are quite complex and depend upon many factors such as geographical localization, environmental factors, management/agricultural practices, and natural enemies that as a whole may be referred as the selection pressures (Boush and Matsumura 1967; Kikuchi et al. 2012; Werren 2012; Kakumanu et al. 2016; Cheng et al. 2017; Gressel 2018; Xia et al. 2018; Vijayakumar et al. 2018).

In addition, it is worth noting that our metagenomic analysis also showed that an important proportion of the reads was assigned to the obligate symbiont “Candidatus Portiera” (a quarter to more than a half) and for the facultative symbionts Hamiltonella and Rickettsia (from one-thousandth to more than two-fifths), for which both “Candidatus Portiera” and Hamiltonella genera belong to Gammaproteobacteria class, and Rickettsia genus belongs to Alphaproteobacteria class; and in turn, these genera become the main constituents of the bacterial community structures from whiteflies B. tabaci MEAM1 natural populations. The importance of the bacterial interaction with insects such as symbionts (obligate and/or facultative) is essential for host survival, since those bacteria are required for the synthesis of essential amino acids and hence allows the insect to feed on a phloem-based diet (Su et al. 2014; Upadhyay et al. 2015; Chen et al. 2016; Liu et al. 2016; Zhang et al. 2016).

Finally, the whiteflies B. tabaci MEAM1 bacterial community structures were not only influenced by their feeding habits of an entirely phloem-based diet, but also the crop plants have an impact over their bacterial symbiont distributions. We were able to determine a slight bias (not significant) on the facultative symbionts (Hamiltonella and Rickettsia) relative abundance in a crop-independent manner (excepting sample 2CH, Solanaceae, Table 3). Interestingly, the obligate symbiont (‘Candidatus Portiera’) practically was not affected by the feeding of different crops of the insect (Table 3). Therefore, the endosymbiont relative abundances from whiteflies B. tabaci MEAM1 natural populations were not affected by the host plant in field conditions, contrasting with previous reports for other whiteflies B. tabaci biotype MED and other aphids at in vitro conditions (Pan et al. 2013; Liu et al. 2016).

Conclusions

The geographic localization and diet of whiteflies natural populations in Baja California Sur had a direct effect to maintain their bacterial community structure, since the main phyla and classes were Proteobacteria, Alphaproteobacteria, and Gammaproteobacteria, respectively, as reported previously for other whiteflies around the world. However, agricultural practices (agrochemicals application) were a factor that influenced on the bacterial community structure at family taxonomic level, and also was reflected in different clustering analyses. The insect host plants also did not have an impact on to the endosymbiont relative abundances, in spite of sap nutrient differences. Finally, to achieve a deeper understanding of the bacterial community structures of whiteflies and the functional role that those carried out in conjunction with the symbionts inside this insect pest, further investigation is necessary.