Molecular and phenotypic variations in Eutetranychus orientalis (Klein) populations from Saudi Arabia

The oriental red spider mite, Eutetranychus orientalis (Klein) is a major pest of citrus in many countries including Saudi Arabia (SA). The morphological variations among the different populations of E. orientalis were reported. In the present study, phenotypic variations based on 40 different morphological characteristics were evaluated in 10 E. orientalis populations collected from different hosts and regions of SA. Further, ITS2-rDNA sequences were used to confirm the identity of these phenotypically varying populations. Phenotypic variations in all populations were found in the shape and length of dorsal setae, striation pattern between setae d1 and e1, and leg chaetotaxy. The rDNA sequence analysis of these morphologically varying populations resulted in 10 different ITS2 Saudi haplotypes. The phenotypic and genetic variations were more related to the host plants rather than their geographic distribution. The E. orientalis population collected from Phoenix dactylifera was phenotypically distinct and genetically divergent. The populations collected from citrus species were also more phenotypically and genetically related to each other than to populations collected from non-citrus host plants. The haplotypes recovered from Ziziphus sp., Morus sp., and Azadirecta indica from different regions were grouped in the same sub-clade. Further, the ITS2 haplotypes of Saudi E. orientalis recovered from Citrus reticulata from Riyadh and Al Ula were 100% identical to the ITS2 haplotypes recovered from Citrus sp. from Israel. It is concluded that phenotypic variations exist among different populations of E. orientalis inhabiting different host plants. This species should be identified carefully by considering phenotypic intraspecific variations.


Introduction
The oriental red spider mite, Eutetranychus orientalis (Klein), is a major pest of many economical shrubs and fruit trees including Citrus sp. and Prunus sp. They are widely distributed in many countries of Oriental, Afrotropical, and Palaearctic regions, including Saudi Arabia (SA) and have been found inhabiting more than 200 host plants [1][2][3][4] variations have been documented among and within different E. orientalis populations [5][6][7][8][9]; these morphological variations in E. orientalis resulted in misidentifications. Recently, the morphological characteristics of the Eutetranychus species were tabulated with descriptions and illustrations; in addition, seven Eutetranychus species as junior synonyms of E. orientalis were suggested based on phenotypic variations [4]. At the intraspecific genetic variation level, the DNA sequences of the internal transcribed spacer 2 (ITS2) were used to genotype nine E. orientalis populations from Israel [10]. It was reported that the morphological and genetic variations in E. orientalis populations could be caused by either vast geographic distribution or polyphagous feeding behavior of the species [11][12][13][14]. Molecular techniques have been used over the past two decades to study the intra / interspecific variations within and among populations / species. Further, DNA-based intra and interspecific variations have been used to identify spider mite species [10,15] and resolve many other cryptic species [16,17]. DNA barcoding using short DNA sequences such as ITS2 or mitochondrial Cytochrome Oxidase I (COI) regions is an effective tool to study intraspecific variations within and among different populations of a species [18][19][20][21].
The aims of the present study were to (1) assess the phenotypic variations in different populations of E. orientalis collected from different localities and host plants in SA and (2) confirm the intraspecific nature of these phenotypic variations through ITS2-rDNA.

Mite collection
This study did not involve any protected or endangered species. Populations of E. orientalis were collected from diverse host plants from different localities in six regions (Jizan, Makkah, Madina, Najran, Riyadh, and Tabuk) of SA during 2009-2018 ( Table 1). The collection sites included agricultural farms and no specific permission was required for the collection of E. orientalis. A distribution map was generated using GPS coordinates of collection sites through ArcGis 10.5, Esri.com computer software (Fig 1). These mite specimens were collected by shaking vegetative plant parts over a white piece of paper. The mites visible under a hand-held magnifying lens and moving on the paper were picked using a camel hair brush and preserved in small vials containing 70% ethanol. These samples were brought to the Acarology laboratory, Department of Plant Protection, College of Food and Agriculture Sciences, King Saud University. Further, few adult females from each population were separated and preserved in 96% ethanol for molecular identification. Also, some adult mite females from each population were mounted on glass slides in Hoyer's medium under a stereomicroscope (SZX10, Olympus, Tokyo, Japan). The species was identified under a phase contrast microscope (BX51, Olym-pus1, Japan) using the identification key [4].

DNA extraction and amplification
Total genomic DNA was extracted from a single female, representing a population collected from 10 different locations (including seven different host plants), using a DNeasy mini kit (Qiagen Venlo, Netherlands). The data of the collection localities and hosts are provided in Table 1.Each extracted DNA sample was properly labelled with the name of the locality and host plant and the date of collection, and directly stored at -20˚C. The concentration of total genomic DNA was evaluated using a Nanodrop spectrophotometer (Thermofisher Scientific, Waltham, Massachusetts, USA) and adjusted for PCR amplification. The PCR primers ITS2forward (5'-GTCACATCTGTCTGAGAGTTGAGA-3') and ITS2-reverse (5'-GTARCCTC ACCTRMTCTGAGATC-3') were used to amplify the ITS2 region of rDNA (Ben David et al. 2007).
The PCR reaction was performed in a 30 μL reaction volume with 15 μL GoTaq 1 Green Master Mix (Promega, Madison, Wisconsin, USA), 0.75 μL of each 10 μM primer, 2 μL (approximately 20 ng) DNA template and 11.5 μL Nuclease free water (Promega). The PCR cycling conditions were as follows: denaturation step for 5 min at 94˚C, followed by 35 cycles of 60 s at 94˚C, annealing at 52˚C for 90 s and 72˚C for 60 s, and the final extension step at 72˚C for 10 min. The PCR products were assayed in 1.2% agarose gel stained with Acridine Orange dye in 1× TAE buffer. These gels were later visualized under UV light using gel documentation system. Further, the PCR products were purified using a MoleQule-On PCR product purification kit (MoleQule On, Auckland, New Zealand).

DNA sequencing and analysis
The purified PCR products were directly sequenced using the same primers used in the PCR amplification; the sequences were then run on ABI 3500 (Life Technologies, MD, USA). The obtained sequences were cleaned and analyzed using the BioEdit software [22]. Thereafter, the cleaned sequences were searched in the NCBI GenBank Database using BLAST. The ITS2 sequences of E. orientalis obtained in the present study and those retrieved from GenBank were aligned using a Clustal W multiple alignment tool. The final sequences were submitted to the NCBI database; the accession numbers of each population are provided in Table 1.
Phenotypic and phylogenetic analyses. The morphological variations between the collected populations were assessed based on different characteristics including length and shape of setae, length and chaetotaxy of legs, and dorsal striation patterns and shape of stylophore. Based on these phenotypic character, an UPGMA dendrogram was constructed using a Minitab 18 data analysis software (2010) (State College, PA: Minitab, Inc. (www.minitab.com). To assess genetic variations among 10 Saudi E. orientalis ITS2 sequences, neighbor-joining (NJ) trees were constructed using the MEGA X tool. The ITS2 haplotypes were determined manually by aligning the DNA sequences using BioEdit software (http://www.mbio.ncsu.edu/ BioEdit/bioedit.html). Single nucleotide polymorphisms (SNPs) were detected and each ITS2 haplotype had its unique SNPs. Similarly, another NJ tree was constructed to assess the genetic variability among the Saudi haplotypes and 12 ITS2 E. orientalis sequences retrieved from the NCBI database. The evolutionary distances were computed using the Tamura-Nei method [23]. To support NJ tree topology, bootstrap analysis [24] with 1000 replications was conducted. The ambiguous positions for each sequence pair were removed through the pairwise deletion method.

Results
Most variations observed in E. orientalis, collected from seven different host plants in 10 locations representing 5 regions of SA, were in the lengths and shapes of dorsal body setae especially dorso-central setae, striation pattern between setae d1 and e1, and the length of legs and the chaetotaxy of leg segments; femora and tibiae (Tables 2 and 3).   In the present study, a UPGMA dendrogram generated from 40 different morphological characteristics (Fig 2A) and an NJ tree constructed based on 10 nucleotide sequences (Fig 2B) divided 10 E. orientalis populations into two main clades and nine sub clades based on host plants (monocot or dicot) rather than their geographic distribution. In both the UPGMA and NJ phylogenetic trees clade 1 contained one E. orientalis phenotype/haplotype, found on Phoenix dactylifera (Arecaceae) (a monocot) from Riyadh, with lowest phenotypic similarity percentage (6.26) and genetically highly distant from other nine haplotypes (Fig 2A and 2B). This phenotype was morphologically distinct from others by having dorso-central setae were comparatively longer and slender in shape ( Table 3).
The clade 2 comprised nine phenotypes/haplotypes of E. orientalis collected from only dicot plants. Among them, four were collected from non-citrus host plants, which were clustered together in both UPGMA (sharing 10 phenotypic characters) and NJ (at 67 bootstrap value) phylogenetic trees (Fig 2A and 2B, Table 3). The phenotypes P4 (from Azadirachta indica (Meliaceae) in Qaseem) and P9 (from Citrus sp. (Rutaceae) grown within A. indica and Morus sp. (Moraceae) trees in Domat ul Jandal) showed highest affinity forming a cluster at 31.67% similarity level (Fig 2A). These two phenotypes shared half the morphological characteristics, which were used to study the phenotypic variations in the present study (Table 2). Similarly, phenotypes P5 and P6 found on Ziziphus sp. (Rhamnaceae) and Morus sp. from

PLOS ONE
Wadi e Dawasir and Qaseem regions, respectively, joined a cluster at 22.42% phenotypic similarity level (Fig 2A) while sharing 16 phenotypic characteristics (Table 3). Among the other five E. orientalis phenotypes/haplotypes collected from citrus host, the phenotypes P1, P8, and P10 joined in a separate sub-clade (sharing only six phenotypic characteristics), whereas P7 and P3 clustered together (sharing 13 phenotypic characteristics) showing phenotypic similarity percentage of 29.05% in the UPGMA dendrogram (Fig 2A, Table 3). However, genetically these were clustered together at 73% bootstrap value except haplotype H8, which was distantly related to all other dicot haplotypes in clade 2 ( Fig 2B).
The ITS2 sequence analysis revealed that the 10 E. orientalis phenotypes gave 10 different ITS2 haplotypes. The length of ITS2 sequence ranged from 496 to 540 bps (Table 1) with an average of 38.8% G-C content. In addition, nucleotide homology among 10 haplotypes was Table 3. Matrix data of morphological character states for Saudi E. orientalis populations. 90-99%. The intraspecific divergence among these haplotypes ranged from 0.0020 to 0.025 (Table 4). The haplotype H3 recovered from C. reticulata in Riyadh was the least divergent (0.0020) from both H1 recovered from Ficus sp. planted in citrus orchard from Al-Ula (Fig 2B,  Table 4) and H10 recovered from C. reticulata from Tabuk. The NJ phylogenetic tree based on 22 ITS2 sequences of E. orientalis haplotypes from different countries including those recovered in the present study showed that most haplotypes were generally clustered based on host plant association rather than geographic distribution. (Fig 3). Haplotypes H1 and H3 were 100% identical to E. orientalis ITS2 haplotypes DQ656471 and DQ656474, previously collected from Israel on Citrus sp., respectively (Table 4, Fig 3). Further, two E. orientalis haplotypes, HQ688670 from Iran (host plant unknown) and DQ656473 from Israel, were 100% identical to each other (Table 4, Fig 3). The haplotype H2 from P. dactylifera (Arecaceae), which was genetically divergent (0.016-0.025) from all other E. orientalis haplotypes, clustered with a Spanish haplotype (KP642048) from citrus agro-ecosystem in the NJ phylogenetic tree (Table 4, Fig 3).

Discussion
In the present study, the phenotypic variations among 10 Saudi E. orientalis populations collected from different hosts and regions of SA were found mostly in lengths and shapes of dorsal body setae especially dorso-central setae, striation pattern between setae d1 and e1, and length of legs and chaetotaxy of leg segments; femora and tibiae (Tables 2 and 3). These variations have been well documented in different E. orientalis populations from different countries [5,10] and have resulted in the addition of species in the genus Eutetranychus [4]. Because of these intraspecific phenotypic variations in E. orientalis populations, four Eutetranychus species viz. E. ricini [25], E. monodi [26], E. sudanicus [27], and E. annecki [1,7,8] were synonymized with E. orientalis. Recently, a study was conducted [4] on the phenotypic variations among different populations of E. orientalis from SA and it reported seven species, viz. E. phaseoli [28], E. guangdongensis [29], E. xianensis [29], E. fici [8], E. pruni [30], E. ricinus [30], and E. sanaae [30] as junior synonyms of E. orientalis. The present study fortifies the findings of a previous study [4] and the comprehensive diagnostic key provided should be followed for the proper identification of Eutetranychus species.
The present study revealed that different hosts associated with E. orientalis races could be present, which might cause phenotypic and genetic differences in these races. Further, E. orientalis is a polyphagous species reported from > 223 diversified host plants including shrubs and trees from both mono and dicot plant families [1]. Host-associated phenotypes have frequently been found in other phytophagous arthropods mainly pest insects [31][32][33] and spider mites i.e. Tetranychus urticae Koch [34][35][36]13] and Panonychus citri (McGregor) [37]. Populations from closely related hosts show closer relationships. It has been reported that the species of spider mites, like T. urticae, T. kanzawai, and Oligonychus gotohi adapted to different hosts were reproductively isolated and genetically different [35,38,39]. Similar results were reported in a previous study [36], that T. urticae populations from the same host Viola sp. from Italy and France were more closely related when compared with those from different hosts from the same country. The species, Brevipalpus phoenicis, which is a host specific species, failed to adapt on alternate hosts and did not exhibit host based genetic variations [40]. However, this species had significant genetic variability among/between populations when collected and compared from different localities in Mexico and Brazil [41].
The utilization of different host plants may cause sympatric speciation in arthropods [42]. The phytophagous species and plants exert strong selective pressures on each other [43]. In addition, herbivores when exposed to and colonized on a new host species, resist and/or Table 4. Genetic divergence (p-distance %) of world E. orientalis ITS2 sequences reported in present study and from NCBI database.  tolerate defense mechanisms specific to that new host, which results in either adaptation or local extinction of the phytophagous species due to strong selection pressure [44]. Spider mites inhabiting wide host range like T. urticae may form host races to combat toxic chemicals produced by different hosts influencing high selection pressure on herbivores to provide opportunity for prolonged existence on different hosts and strong association [45].
In the present study, the genetic data of 10 Saudi E. orientalis populations was associated with host plants they inhabited (Fig 2A and 2B). The phylogenetic tree based on world E. orientalis ITS2 sequences, including from the present study, supported host association (Fig 3). Similar findings were reported in Tetranychus urticae (Koch) [36] and Amphitetranychus viennensis (Zacher) [46]. In a previous study, nine E. orientalis populations from Israel have been reported as genetically distinct [10]. The intraspecific genetic divergence of nine Saudi E. orientalis haplotypes was in a similar range (0.002 to 0.017) as found in the Israeli E. orientalis haplotypes [10].
The Saudi haplotype H1 collected from Ficus sp. in a citrus orchard was 100% identical to Israeli haplotype (DQ656474 from Citrus sp.) and shared the clade with other Saudi haplotypes collected from Citrus sp. Likewise, the genotype H9 from Domat ul Jandal on Citrus sp. under the A. indica trees was clustered in the phylogenetic tree and dendrogram with G4 from A. indica from Qaseem; this might be due to the colonization of E. orientalis on a new host. The

PLOS ONE
herbivores adapt on new hosts with the passage of time after coping with the defense mechanisms of new hosts [44,47]. The wide range of host plants and geographic distribution of E. orientalis indicate that this species has a strong ability to adopt to diversified hosts and cope with defenses exerted by the host species, thus making it possible to live successfully in different environmental conditions and topography. The mite species occurring on different host plants undergo marked changes in their behavior, morphology, and genetic structure. This commonly happens in polyphagous mite species where adaptation to different host plants results in genetic differentiation [48] In addition, the intraspecific genetic variations in mites could also be due to geographic distances. The geographic distances act as a barrier for gene flow and promote intraspecific differentiations. The mite species, Lepidoglyphus destructor, which have active dispersal capabilities, have shown genetic differentiations among populations collected from localities that are far apart from each other [49]. In the case of grapevine bud mites, which have limited dispersal, strong genetic differentiation was reported among the specimen collected from different vineyards and even different parts of same vineyard [50]. The genetic diversity was also recorded in the Citrus mite, Panonychus citri, collected from different localities in China [51].

Conclusion
The polyphagous behavior and wide distribution of E. orientalis resulted in variations in phenotypic characteristics as well as genetic structures. The taxonomists should consider the phenotypic variations specially in the shape and length of dorsal setae and striation pattern between dorsal setae d1 and e1, while identifying different E. orientalis populations to avoid improper identifications. Further extensive studies are needed to investigate the occurrence of different races of E. orientalis associated with a specific host or locality. Also, its population structure needs consideration in different agro-ecosystems.