Developing conversed microsatellite markers and their implications in evolutionary analysis of the Bemisia tabaci complex

The study of population genetics among the Bemisia tabaci complex is limited due to the lack of conserved molecular markers. In this study, 358, 433 and 322 new polynucleotide microsatellites are separately identified from the transcriptome sequences of three cryptic species of the B. tabaci complex. The cross species transferability of 57 microsatellites was then experimentally validated. The results indicate that these markers are conserved and have high inter-taxon transferability. Thirteen markers were employed to assess the genetic relationships among six cryptic species of the B. tabaci complex. To our surprise, the inferred phylogeny was consistent with that of mitochondrial COI sequences, indicating that microsatellites have the potential to distinguish species of the B. tabaci complex. Our results demonstrate that development of microsatellites from transcriptome data is a fast and cost-effective approach. These markers can be used to analyze the population genetics and evolutionary patterns of the B. tabaci complex.

To date, 54 microsatellite markers are available for B. tabaci 14,[21][22][23] . However, all of these microsatellites were derived from the genomic DNA and the connections between these markers and gene functions are completely unknown. Expressed Sequence Tag (EST) and transcriptome sequences contain polymorphic genetic markers and can be used to identify microsatellites 24 . Compared to the genomic DNAderived microsatellites, EST-and transcriptome-derived microsatellites lack introns and intragenic regions and can correspond to genes with known or predicted functions 25 . In addition, those microsatellites have fewer null alleles and stutter bands 26 and have more potential statistical power in multiple comparisons 27 . Furthermore, EST-and transcriptome-derived microsatellites have high degree of transferability across species 28 , and can be used in closely relative species 29 . Therefore, systematical investigation of EST-and transcriptomederived microsatellites will facilitate evolutionary and comparative studies in the B. tabaci complex composed of closely related cryptic species 27 .
Recently, the transcriptomes of two invasive (MEAM1, MED) and one indigenous whitefly species (Asia II 3) have been sequenced 30,31 . These studies have generated a tremendous amount of data and provided a valuable source for the identification of microsatellite markers in whiteflies. The first objective of this study is to identify microsatellites from the three transcriptomes. In addition, microsatellites located in different regions of a gene serve various functions 32 . The distribution of microsatellites on genes was also analyzed. Furthermore, PCR experiments were employed to verify these predicted microsatellites and their cross species transferability. By comparative analysis of the newly developed microsatellites, the genetic relationships of six B. tabaci species were revealed. This study provides a rich resource of microsatellites for the B. tabaci complex and will facilitate researches on whitefly genetic diversity and evolution.

Results
Identification of microsatellites from the B. tabaci transcriptome databases. A total of 27.653 Mbp, 44.937 Mbp and 24.468 Mbp of sequences from the MEAM1, MED and Asia II 3 transcriptomes were used for mining microsatellites with the MISA-Micro Satellite program 33 (Table 1). There were 6419, 11711 and 4115 microsatellites in MEAM1, MED and Asia II 3 respectively (Table  S1), which correspond to one microsatellite per 3.837 , 5.946 Kbp of transcriptome sequences. The total numbers of polynucleotide repeats were 358, 433 and 322 in MEAM1, MED and Asia II 3 respectively. While most microsatellites-containing unigenes have only one microsatellite, there are 362, 299 and 190 unigenes containing multiple microsatellites. In addition, 277, 367 and 193 unigenes contain compound microsatellites in MEAM1, MED and Asia II 3, respectively (Table 1).
Of the characterized microsatellites, mononucleotide repeats were the most common, followed by dinucleotide, trinucleotide and tetranucleotide repeats (Fig. 1A). On a complementary strand, a polyA repeat is the same as a polyT repeat. Similarly, in different reading frames or on a complementary strand, (AC)n is the same as (CA)n, (TG)n and (GT)n, while (AAG)n is the same as (AGA)n, (GAA)n, (CTT)n, (TTC)n, and (TCT)n. Thus, mononucleotide, dinucleotide and trinucleotide repeats can be grouped into 2, 4 and 10 unique classes respectively 34 . In the three species of the B. tabaci complex, A/ T motifs were the most abundant in mononucleotide repeats, and AG class were the most common in dinucleotide repeats (Fig. 1B). However, the usage of trinucleotide was different. In MEAM1, the AAG class was the most widespread followed by ATG and AAC class (Fig. 1C). In MED, the most prevalent three triplet codons were AAG class, whereas in Asia II 3, ATG is the most frequent class (Fig. 1C). A total of 45 tetra microsatellites were identified from both MEAM1 and Asia II 3. Furthermore, two pentra motifs were found in Asia II 3 and one hexa motif was found in MEAM1 and Asia II 3, respectively (Fig. 1A).
Distribution of microsatellites in 39UTR, 59UTR and CDS regions. The distribution of polynucleotide microsatellites in CDS, 59UTR and 39UTR regions was investigated. Based on the information of BLASTx homology, the position of 71, 71 and 40 polynucleotide microsatellites were respectively determined for MEAM1, MED and Asia II 3 ( Table 2 & Table S2). In the CDS, 45, 48 and 27 microsatellites with polynucleotide repeats were found in MEAM1, MED and Asia II 3, respectively, which were significantly higher than that of UTRs (Table 2). Interestingly, the number of trinucleotide repeats in CDS region was also much higher than other types of microsatellites. The characteristics of the amino acids encoded by the trinucleotide repeats in CDS were then investigated. In MEAM1, MED and Asia II 3, a total of 37, 37 and 20 triplet codons were found and they encoded 10, 15 and 8 different amino acids respectively ( Fig. 2A). The codons encoded aromatic amino acids took up most of the partition, followed by aliphatic and heterocyclic amino acids (Fig. 2B). The codons encoded hydrophilic amino acids were 29, 20 and 17, while encoded hydrophobic amino acids were 3, 15 and 3 in MEAM1, MED and Asia II 3, respectively.
Gene Ontology (GO) and KEGG annotation of microsatellitecontaining sequences. GO assignments were used to classify the functions of the genes with microsatellites. Based on sequence homology, 272, 456 and 255 microsatellite-containing sequences from MEAM1, MED and Asia II 3 have GO annotations and can be categorized into 35 functional groups. The 35 functional groups were classified into three main categories (Biological process, Cellular component and Molecular function) (Fig. 3). GO analysis showed that 'Metabolic process', 'Cellular process' and 'Cell part' terms are dominant. Next, these genes were annotated to different KEGG pathways. A total of 318, 538, 297 microsatellite-containing sequences of MEAM1, MED and Asia II 3 were mapped to 193, 222 and 210 KEGG pathways, respectively (Table S3). Some of these genes are related to resistance to environmental stresses and insecticides, such as aldehyde oxidase 35 , cytochrome P450 36 and mitogen-stress activated protein kinases 2 37 (Table S4).
Characterization of the predicted microsatellite markers. To validate the predicted microsatellites, 88 primer pairs were synthesized (24 for MEAM1, 32 for MED and 32 for Asia II 3) to amplify microsatellites from whitefly DNA. Among the 88 primer pairs, 57 (65%) generated clear PCR products and 31 (35%) did not  Table S5.
The cross-taxa transferability of microsatellites. Next, the crossamplification of the 57 primer pairs were investigated by PCR and capillary electrophoresis 38 . The result showed that 42 (73.68%) amplified fragments from all of the three whitefly species (Table  S6), suggesting that these microsatellites are highly conserved and may act as genetic markers for the B. tabaci complex. In addition, 8 primer pairs could amplify fragments from the two invasive whiteflies and 1 pair amplified from both MED and Asia II 3, while the remaining 6 pairs of primers only amplified fragments from one of the three species. Next, the functions of these 42 microsatellite-containing genes were determined through Blast search, of which 33 showed significant similarities to known genes (Table S5). The alleles of these 42 microsatellites were compared among MEAM1, MED and Asia II 3 whiteflies. Of the 42 microsatellites, 26 contained the same alleles in MEAM1 and MED, whereas only 17 markers shared in MEAM1 and Asia II 3 and 16 markers shared in MED and Asia II 3 (Table S6), which is consistent with the fact that MEAM1 and MED have closer relationships among the three cryptic species. Interestingly, the sodium channel gene (Gene ID 99267), which is associated with xenobiotic resistance 39 , was found to contain different alleles in invasive and native whitefly species (Table S6).
Characteristics of the 13 microsatellites in six species of the B. tabaci complex. Thirteen polymorphic microsatellites (Table S5 marked in yellow) were then employed to assess the polymorphism and heterozygosity among six laboratory colonies of the B. tabaci complex (8 individuals for each species) (    glycine-rich cell wall structural protein precursor (Gene ID: 29116), sodium channel (Gene ID: 99267) and NADH dehydrogenase subunit 4L (Gene ID: 76476) had high polymorphism ( Table 3). The CCAAT/ enhancer binding protein (Gene ID: 102573) showed relatively low polymorphism ( Table 3). The polymorphism information content (PIC) ranged from 0.246 to 0.866, with an average of 0.612 in the six cryptic species, which indicates the effectiveness of microsatellites markers for detecting polymorphism 40 . In addition, the genetic diversity of the 13 microsatellite loci was compared among the six whitefly species. The level of gene diversity from the highest to the lowest was Asia II 3, MED, China 1, Asia II 6, MEAM1 and Asia II 7 (Table S7). Fig. 4 displays the gene diversity of the 13 microsatellite loci in every comparison. For example, in the comparison of MED and MEAM1, 7 of the 13 loci in MEAM1 showed a decrease of gene diversity compared to that of MED. To further testify these markers, 4 microsatellites were randomly picked to sequence and analyze the alleles per locus among the six cryptic species (Fig. 5). Complex mutational patterns (single base mutation, change of repeat units and indels within flanking region) were observed in these transcriptomederived microsatellites (Fig. 5), which is normal in insects 41 . . When considered the first and third factors, the results revealed clearly that there are three genetic groups among the 6 cryptic species. The PCA analysis supports the phylogenetic clusters. In conclusion, these microsatellites may be used as markers to describe the genetic diversity of the B. tabaci complex.

Discussion
To our knowledge, this is the first investigation of the frequency and distribution of B. tabaci microsatellites using transcriptome data (Table 1). We found that mononucleotide repeats were the most common, followed by dinucleotide repeats, trinucleotide and tetra-   nucleotide repeats in all of the three species. These results confirm the theory that the microsatellite abundance decreases with the increase of motif repeat number and repeat length 43 . These microsatellites provide a valuable resource for the development of genetic markers in B. tabaci. In addition, these markers offer a chance to classify the functions of these microsatellite-containing genes. Interestingly, some microsatellite-containing genes were found in pathways related to environmental stress responses (Table S4).
These markers may open a new avenue for the research on the B. tabaci pesticide resistance by estimating the frequencies of alleles in genes related to resistance across populations 44 .
Many studies have demonstrated that microsatellites derived from transcribed sequences harbor higher transferability because of their conservative features 45,46 . In this study, among the 57 microsatellite markers, 42 primer pairs (73.68%) could amplify fragments from the three species. These results illustrate that the transcriptome-derived markers have high inter taxon transferability 47 . Therefore, these markers can be used to amplify microsatellites from the other closely related species that do not yet have markers and can be considered as anchor markers for comparative mapping and evolutionary studies across species 48 .  ity are critical characteristics of microsatellite markers. The high allele numbers and heterozygosity (Table 3) of these transcriptome-derived markers suggests that they can be used to assess the genetic profile of the B. tabaci complex. A total of 25 loci significantly deviated from Hardy-Weinberg equilibrium ( Table 3). The source of deficiencies includes the recurrent inbreeding, subpopulation structure (Wahlund effect) and/or null alleles. The big differences between H o and H E also were found among the genomics-based microsatellites developed for B. tabaci by De Barro 49 . In addition, the marked differences between H o and H E also were observed in the genomicsbased microsatellites of B. tabaci collected in different environment 14,16 . Some loci did not have heterozygosity and the F IS values were all positive, suggesting that the loss of heterozygosity can happen in some microsatellites. This is probably due to the fact that our samples were obtained from lab populations. Additional experiments with field populations are warranted to reveal the polymorphism and heterozygosity of different species of the B. tabaci complex. Due to the lack of morphological characteristics, the systematics of B. tabaci species complex depends largely on molecular approaches [2][3][4]6 . To date, the evolutionary history of B. tabaci complex was inferred exclusively from the partial sequence of COI gene 2,6,50 . Other molecular markers are in great need to further illustrate the complicated evolutionary relationships of the B. tabaci species complex. For the 4 major clusters that represent the B. tabaci species complex, the most divergent cluster is SubSaharan Africa followed by New World. We have used 11 newly developed microsatellites to study the SubSaharan Africa1 and the results revealed that they can work on the cryptic species (data unpublished). Owing to the absence of the specimen, we did not study the New World, however, the most divergent cluster is SubSaharan Africa followed by New World among the 4 major clusters. Therefore, our newly developed microsatellites are conserved and many of them can be used to analyze other cryptic species of the B. tabaci complex. The genetic relationship of 6 species derived by 13 novel conserved microsatellites markers is in accordance with the phylogeny of the partial COI sequences. It provides additional evidence for the evolutionary relationship of the 6 cryptic species. In addition, our results also prove that microsatellites can be considered as suitable markers for future evolutionary analysis in other species complex.

Conclusion
In this study, we characterize a large number of transcriptomederived microsatellites from members of the B. tabaci complex. Functional categorization of these microsatellites-containing genes provides valuable information about the potential functions of these microsatellites. In addition, 57 microsatellite markers were experimentally validated. Many of these microsatellites can be used for population analysis and are transferable across different species of the B. tabaci complex. Complex mutational patterns were observed in these transcriptome-derived microsatellites. What's more, by ana-  lyzing the genetic relationships of six B. tabaci species basing on newly developed microsatellites, the results indicated that these microsatellites may be used as markers to describe the genetic diversity of the B. tabaci complex. These markers enrich the existing microsatellite markers of B. tabaci and can be used to analyze the genetic diversity and evolutionary pattern of this whitefly complex.

Methods
Mining the microsatellites from the transcriptome database. The MISA-Micro Satellite identification tool (http://pgrc.ipk-gatersleben.de/misa/) was used to search for microsatellites in the transcriptome data of MEAM1, MED and Asia II 3 whiteflies. Microsatellites were defined as being mononucleotide repeats .5 10 repeats and di-, tri-, tetra-, penta-and hexanucleotide repeats .5 6 repeats 33 . Criteria for compound microsatellites was an interval of bases ,5 100 of the motif length.
The location of microsatellites. Coding sequences (CDS) of each gene were first determined by BLASTx against the Swissprot database using a threshold of 1 3 10 25 .
CDS with unexpected stop codon in the Blast hit region were removed. Start codon positions were determined by examination of the in-frame ATG codon present 30 bp upstream or downstream of the beginning of the aligned reference protein.