Molecular mechanisms of Tetranychus urticae chemical adaptation in hop fields

The two-spotted spider mite, Tetranychus urticae Koch is a major pest that feeds on >1,100 plant species. Many perennial crops including hop (Humulus lupulus) are routinely plagued by T. urticae infestations. Hop is a specialty crop in Pacific Northwest states, where 99% of all U.S. hops are produced. To suppress T. urticae, growers often apply various acaricides. Unfortunately T. urticae has been documented to quickly develop resistance to these acaricides which directly cause control failures. Here, we investigated resistance ratios and distribution of multiple resistance-associated mutations in field collected T. urticae samples compared with a susceptible population. Our research revealed that a mutation in the cytochrome b gene (G126S) in 35% tested T. urticae populations and a mutation in the voltage-gated sodium channel gene (F1538I) in 66.7% populations may contribute resistance to bifenazate and bifenthrin, respectively. No mutations were detected in Glutamate-gated chloride channel subunits tested, suggesting target site insensitivity may not be important in our hop T. urticae resistance to abamectin. However, P450-mediated detoxification was observed and is a putative mechanism for abamectin resistance. Molecular mechanisms of T. urticae chemical adaptation in hopyards is imperative new information that will help growers develop effective and sustainable management strategies.

. The acaricide spray model at hopyards during hop season in 2013. Several acaricides with different mode of actions were applied to control T. urticae. Among them, abamectin, bifenazate, and bifenthrin were commercially important acaricides used in hopyards 10 .
Scientific RepoRts | 5:17090 | DOI: 10.1038/srep17090 Acaricide resistance levels in field populations. The toxicities of abamectin and bifenazate were assessed for T. urticae populations collected from 13 and 12 hopyards, respectively (Tables 2 and 3). In the bioassays with abamectin, the LC 50 s ranged from 1.36 to 26.05 mg a.i./L and the resistant ratios (RRs) compared with the susceptible strain varied from 5.96 to 114.25 (Table 2). Low resistance levels (RR < 10) were observed in 10.5% of the surveyed populations, 10.5% had high resistance (RR > 100), and the majority of the surveyed populations (79%) exhibited moderate resistance (RR = 10-100) to abamectin (Fig. 3A). The RR of the T. urticae population in the organic hopyard (Granger 2) compared with the susceptible population was 11.23, which is the 3 rd lowest resistance among surveyed populations and the highest level of mortality (100%) at the field rate. Samples collected from the Granger 4 hopyard showed the lowest resistance ratio (RR = 5.96) compared with the susceptible population. There were three 1 st year (baby) hopyards (Prosser 3, 4, and 5) surveyed in 2013. The RRs of samples collected from these baby hopyards ranged from 21.80 to 114.25, exhibiting a moderate to high degree of resistance (Table 2). There were multiple collections from certain hopyards (Prosser 2, 3 and 4) during the course of summer 2013. Specifically, six collections were taken from the Prosser 2 hopyard starting from middle of June till just prior to harvest in late August during which abamectin was applied twice 10 . The RR increased 6-fold from the middle of July to mid-August ( Table 2). The RRs in samples collected from Prosser 3 and 4 increased 1.7-fold and 2.3-fold in four and five weeks, respectively. The highest resistance level to abamectin was recorded at the Prosser 4 (RR = 114.25) ( Table 2).
In the bioassays with bifenazate, the LC 50 s ranged from 3.93 to 78.97 mg a.i./L and the RRs varied from 4.79 to 96.30 (Table 3). Populations exhibiting low resistance levels (RR < 10) accounted for 37.5% of the populations surveyed, and 62.5% of the populations exhibited moderate resistance (RR = 10-100) to bifenazate (Fig. 3B). The lowest RR to bifenazate, 4.79, was recorded from the samples collected from the organic hopyard (Granger 2). The RRs of samples collected from the 1 st year hopyards showed low to moderate level of resistance ( Table 3). The highest RR to bifenazate was observed in the sample collected from Granger 3 (RR = 96.30) ( Table 3). Due to the limited number of collected T. urticae individuals in four populations, only the discriminating dose of bifenazate was evaluated ( Table 3).

Evaluation of target site mutations. The occurrence of 16 mutations in four target genes, GluCl1
and GluCl3 (target of abamectin; Fig. S1), cytb (target of bifenazate; Fig. S2), and VGSC (target of bifenthrin; Fig. S3), was examined in T. urticae field populations by direct sequencing of PCR products. By visual examination of sequencing chromatographs at the mutation sites, we could identify samples that contained wild-type, resistant, or both alleles. The combination of mutations in field T. urticae populations collected from PNW hopyards exhibited a unique pattern (Table 4). Only two mutations, G126S  and F1538I, in cytb and domain III of VGSC, respectively, were identified (Table 4). There were no mutations observed in GluCl1, GluCl3, and other region of cytb and VGSC.
No mutations observed in Glutamate-gated chloride channel genes. Inhibitory Glutamate-gated chloride channels (GluCls), members of the cys-loop ligand-gated ion channel (cysLGIC) superfamily, are extrajunctional or postsynaptic receptors found in muscle or neural ganglion of most protostome phyla including Chelicerates such as T. urticae 23,33 . The genome of T. urticae contains six orthologous GluCl genes 23 . Previous studies revealed that two mutations in two different GluCl channel subunits, GluCl1 and GluCl3, were related to abamectin resistance in T. urticae 23,26 . Thus we designed primers to sequence the fragments containing these two mutations ( Fig. S1) from susceptible and all hop field populations of T. urticae. Surprisingly, there were no mutations identified from the samples tested (Table 4), suggesting target site insensitivity-mediated resistance is not the mechanism leading to the abamectin resistance that we observed in T. urticae field populations.

Identification of mutations in the cytb gene.
Recent studies suggested that bifenazate resistance was closely correlated with mutation(s) in the mitochondrial cytb 27 . A combination of at least two cd1 helix mutations in the Q o pocket (G126S and I136T or G126S and S141F) and one mutation in the ef helix of Q o pocket (P262T) were linked with a high level of bifenazate resistance in T. urticae. We sequenced an 828 bp fragment of the T. urticae cytb gene, which included the G126, I136, S141, D161 and P262 sites ( Fig. S2) that have been demonstrated to confer bifenazate resistance in T. urticae 27 . One amino acid substitution, G126S, was detected in T. urticae field populations. 35% of field samples field samples contained only the resistant allele, 20% contained both alleles (G/S) and 15% only the susceptible allele (G) ( Table 4; Fig. 4A). Since the G126S mutation alone only causes low to moderate bifenazate resistance 27 , this result is consistent with the bifenazate resistance phenotype observed ( Table 3).
Identification of mutations in the voltage-gated sodium channel gene. The voltage-gated sodium channel (VGSC) is an integral transmembrane protein that is responsible for the rapidly rising phase of action potentials on the neuronal membranes. Due to its essential role in electrical signaling, VGSC is the target of several neurotoxins, including pyrethroids and DDT 34 . Many amino acid substitutions associated with pyrethroid resistance in arthropods are located in transmembrane segments 4-6 of domain II (IIS4-IIS6) including M918 (super kdr), L925, T929, L932, V1010, L1014 (kdr), and L1024 30,34-36 . One mutation within the intracellular inter linker connecting domains II and III (A1215D) and one mutation in domain III (F1538I) were detected in a highly bifenthrin resistant T. urticae strain from Greece 29 . Thus we amplified three fragments of the VGSC from the domain II, II-III inter linker, and domain III regions (Fig. S3). We identified only one amino acid substitution, F1538I. It was observed in 16 out of 24 field samples tested (66.7%), 12 of which contained both alleles (F/I) and 4 of which were only contained the isoleucine substitution (I) ( Table 4; Fig. 4B).
Cytochrome P450-mediated metabolic detoxification. Besides target site insensitivity, cytochrome P450-mediated detoxification had been shown to be one of the most important mechanisms in acaricide resistance of T. urticae [37][38][39] . The genome of T. urticae contains 86 P450 genes. We examined the relative expression of three P450s, CYP385C4, CYP389A1, and CYP392D8, belonging to the CYP3, CYP4, and CYP2 clans, respectively. We chose these three P450s because they have been shown to exhibit more than two-fold up regulation after switching host plants and their expression patterns have been linked to acaricide resistance in T. urticae 9 . The expressions of these three P450s in five field populations from five major locations were compared with their expressions in the susceptible population. As shown in Fig. 5, CYP385C4 had significantly higher expression in all five field populations. However, this increase in expression was not large (less than two-fold). CYP389A1 only showed significantly higher

Discussion
Due to a very short residual effectiveness, abamectin has become the predominant acaricide applied to control T. urticae outbreaks in August as the hops near harvest. Annually, approximately 98% of the hop acreage in Washington is treated with abamectin at least once and 80% is treated at least three times. The widespread use of abamectin on hops raises the distinct possibility of control failure as a result of resistance. From sampling the same hopyard over multiple time points in the same season, we found increasing levels of abamectin resistance, suggesting selection pressure from abamectin applications was driving increasing resistance. For instance, multiple collections in the Prosser 2 hopyard showed that the RR to abamectin increased 6-fold from the middle of July to mid-August ( GluCls together with gamma-aminobutyric acid (GABA)-gated channels and histamine-gated chloride channels (HisCls) are known targets of the macrocyclic lactones, the avermectins (including abamectin) and ivermectins 26,33,43,44 . The point mutation G323D in GluCl1 was tightly linked to a moderate abamectin resistance (17.9-fold) in the AbaR strain 26 . Two point mutations, G323D and G326E, in GluCl1 and GluCl3, respectively, were identified in a > 2,000-fold abamectin resistant strain 23 . However, there was no mutation on GluCl subunits detected in any hop samples (Table 4), suggesting target site insensitivity is not likely the mechanism involved in resistance to abamectin in T. urticae field populations. Our results are comparable with a study by Khajehali et al. 45 which found no GluCl mutations in 15 T. urticae strains collected from rose greenhouses in the Netherlands, although 10 of those strains displayed abamectin resistance. Many recent studies also suggested that abamectin target site mutations are not especially common in T. urticae populations worldwide. For example, the G326E was detected in only seven out of 51 T. urticae populations sampled from 27 countries and five continents 46 . The G323D mutation was only found in two Greek samples in the same survey 46 . In another study with 25 Korean T. urticae populations, only one field-collected T. urticae sample contains G323D mutation 47 . Previous synergism tests and transcriptomic data indicated that additional mechanisms such as enhanced metabolic detoxification by cytochrome P450s may be implicated in the abamectin resistance phenotype 37,48,49 . A genome microarray analysis revealed several cytochrome P450 genes were up-regulated in an abamectin resistant strain 49 . Further evidence confirmed the function of one of these P450s, CYP392A16, in metabolizing abamectin 50 . Unfortunately, this study was published after the completion of our study, and we did not have enough sample material remaining to test for expression of this gene. However, of the three P450s we did examine in our study, one Clan 2 P450, CYP392D8, showed constitutive over-expression in all five field collected samples compared to the susceptible population, indicating its potential function in abamectin resistance (Fig. 5).
Bifenazate is a hydrazine carbazate acaricide that was discovered in 1990 by Uniroyal Chemical and first registered in the state of Washington in 2002 8,51 . Because of the quick knockdown and long residual effects on many economically important phytophagous mite species and low toxicity on predatory mites and beneficial insects, bifenazate is widely used as a selective acaricide to control T. urticae in hopyards. Our bioassay data demonstrated that the majority of field T. urticae populations (62.5%) in hopyards exhibit moderate levels of resistance to bifenazate (Fig. 3B). Our target site mutation screening revealed that a mutation G126S on cytb gene occurs in 35% of T. urticae populations (Fig. 4A). It should be noted that G126S (GGA to AGA) is the same mutation as described in previous studies 27,28 . G126S is the most common substitution on cytb gene of T. urticae that was identified in several bifenazate resistant populations 22,27,28,45,46 . Previous studies showed that mutations on the G137 site in Saccharomyces cerevisiae (equivalent to G126 in T. urticae) contributed to respiratory-deficiency through affecting stability of FeS 52,53 . However, the G126 mutation alone only confers low to moderate level of resistance to bifenazate 28 . In our results, the G126S mutation was observed in populations of Granger 4, Moxee 1 & 2, and Prosser 1, 4 & 5, which all demonstrated low to moderate level of bifenazate resistance (RR = 8.37-23.02 or mortalities at discriminating dose ranged from 90% to 100%) ( Table 3), suggesting the resistance phenotypes of these samples are consistent with their genotypes. Other mutations or mutation combinations on cytb gene such as P262T, G126S with I136T/S141F that are responsible for high bifenazate resistance with RR > 2778 27,28 were not detected in any of our samples.
Bifenthrin, a pyrethroid, has been introduced for T. urticae control in hopyards since it was registered in 1993 51 . Because of their safety, longevity of residual activity and low cost, pyrethroids are extensively used for pest control, with about a 20% insecticide market share 23 . Unfortunately, ubiquitous resistance to pyrethroids had been broadly reported in various insect populations 19,20,34,36,54 . In T. urticae, two mutations, F1538I in domain IIIS6 and A1215D within the intracellular inter linker connecting domains II and III were linked with high bifenthrin resistance in a Greek population 29 . The function of the F1538I mutation in pyrethroid resistance has been confirmed 34,55,56 while the function of A1215D is still unknown. Another substitution, L1024V in domain IIS6 was reported to play an important role in the fenpropathrin resistance of T. urticae from Korea 30 . Pyrethroids are not used very often in the hopyards in PNW since they are linked with subsequent T. urticae outbreaks. Therefore we omitted the toxicity evaluation of field collected T. urticae samples. However, based on our record of acaricide sprays in hopyards, bifenthrin is still used in August as the hops near harvest 8 (Fig. 1). Our DNA diagnostic results demonstrated that a mutation in the VGSC gene (F1538I) was observed in 66.7% T. urticae populations (Fig. 4B). Particularly, F1538I was fixed in four samples collected from Grandview and Prosser 4 & 5, three of which were collected in late August or September (Table 4). Additionally, esterase-mediated metabolic detoxification had also been proposed to confer resistance to pyrethroids in T. urticae [57][58][59] . This result suggests that developing of pyrethroid resistance in hopyards should be of concern.
In summary, T. urticae populations in hopyards exhibit a low to moderate level of acaricide resistance. The mechanisms of acaricide resistance in T. urticae are likely mediated by a number of different pathways: not only target site insensitivities but also enhanced metabolic detoxification. It is a common phenomenon that multiple genes or mechanisms confer resistance simultaneously to a certain pesticide 18,20,[60][61][62][63][64] . Therefore, we plan a genome-wide investigation to identify a more complete set of candidate resistance genes from T. urticae populations of hopyards. Our data also suggests that acaricide spray history, neighboring plants, and time of the season are important factors in correctly diagnosing acaricide resistance in T. urticae. Developing a baseline effective dose for commonly used acaricides and screening local T. urticae populations with resistance-associated molecular markers would be a proactive approach toward T. urticae resistance management. Our study reveals a unique phenotypic and genotypic pattern underpinning the chemical adaptation of T. urticae in hop fields which will be of assistance in developing diagnostic tools for integrated T. urticae management. Mite-infested hop leaves were stored in a plastic bag and transported to the lab in a cooling box within a few hours of collection. Spider mites were identified under a dissecting scope according to morphological characteristics 13 . Approximately 50-100 adults were stored in 95% ethanol for genomic DNA extraction. About 300 adult T. urticae from each of five major locations listed in Table 4 were also stored in RNAlater ® (Sigma-Aldrich, Saint Louis, MO) for RNA extraction. Remaining mites were used for bioassays directly. Because three samples had a low number of mites (Table 4), we reared them on lima bean plants in an isolated walk-in growth chamber for one month to increase population size before sampling them for DNA extraction and bioassays.
Bioassays and data analysis. Leaf disc bioassays were used to estimate the LC 50 (lethal concentration required to kill 50% of the individuals in a population) of abamectin and bifenazate for lab susceptible and field spider mite populations. The method followed that of Knight et al. 65 . Briefly, ten female adult spider mites were placed on the back of a bean leaf disc (2 cm diameter) with a fine brush. Two leaf discs were arranged on water-saturated cotton (4 cm × 4 cm) in a single petri dish (9 cm diameter, 1.5 cm height; Alkali Scientific, Pompano Beach, FL). The water-saturated cotton was pushed up against the perimeter of the leaf disc to prevent mites from walking off the disc 65 . Two commercially formulated acaricides for leaf disc bioassay are Epi-mek ® 0.15 EC (2% a.i. Abamectin, Syngenta Crop Protection) and Acramite ® 50WS (50% a.i. Bifenazate, Chemtura Agro Solutions). The recommended field concentrations for these two acaricides are 23 mg a.i./L and 899 mg a.i./L, respectively. The field rate solutions were prepared in the lab using commercial formulated acaricides and distilled water. These solutions were serially diluted in distilled water for 4-7 concentrations ranged from 0.1-67 mg a.i./L and 0.44-889 mg a.i./L for Epi-mek ® and Acramite ® , respectively.
The sticky tape method was used to estimate the LC 50 to bifenthrin for the lab susceptible strain because pyrethroids are shown to have repellent effects on mites 66 . In this method, ten female adult spider mites were placed dorsal side down on a strip of double-sided sticky Scotch ® tape (3cm × 1.2 cm) stuck on a glass slide (7.5 cm × 2.5 cm). The commercially formulated bifenthrin was Bifenture ® EC, a pyrethroid provided by United Phosphorus (25.1% a.i. Bifenthrin). These bifenthrin solutions were serially diluted in distilled water for 4-7 concentrations ranged from 6-120 mg a.i./L.
Leaf discs or glass slides were treated topically with 2 ml of acaricide solutions with a Potter spray tower (Burkard Manufacturing, Richmansworth, Herts, UK) 67 . The tower was calibrated to deliver 1.1 kg/ cm 2 which allowed 2.0 ± 0.1 mg/cm 2 spray fluid. Each bioassay consisted of 4-7 acaricide concentrations with 4-6 replicates for each concentration. The spider mites exposed to distilled water in the Potter spray tower were used as the non-treated control. The treated leaf discs or glass slides were maintained at 25 ± 2 °C and a photoperiod of 16:8 (L:D) h after the initiation of the bioassay. Mortality was evaluated after 24 h. Mortality was assessed by gently touching each individual spider mite with a fine camel hair paint brush under a dissecting stereomicroscope. The individuals with no response were counted as dead. The few moribund individuals that were not able to maintain balance and show uncoordinated twitching were also recorded as dead. The slope, intercept, and LC 50 (corrected against the untreated control) were evaluated with Abbott's formula 68 calculated by log-dose probit analysis (POLO Probit 2014). The statistical analysis of LC 50 values was based on non-overlapping 95% CI. Resistance ratios (RRs) were calculated through dividing LC 50 values of field samples by the LC 50 value of the lab susceptible population.
Resistance-associated amino acid substitution screening. Genomic DNA was extracted using a DNeasy Blood & Tissue kit (QIAGEN) from 10 adult mites for each population. The DNA was stored at − 20 °C till use. The genomic DNA was used as a template for PCR performed in a Peltier-Effect thermal cycler (MJ Research, Inc., Canada). Primers for PCR amplification of regions with resistance-associated point mutations are listed in Table S1. PCR was performed using Phusion High-Fidelity DNA Polymerase (Thermo Scientific, Pittsburgh, PA) under the following cycling parameters: 95 °C for 3 min 50 s, 35 cycles of 94 °C for 35 s, 55 °C for 35 s, and 72 °C for 3 min, with final extension for 10 min at 72 °C. PCR products were purified using DNA Clean & Concentrator (Zymo Research, Irvine, CA) following the manufacturer's protocol. The purified DNA from each individual was directly sequenced using primers described above (Table S1) for PCR amplification. Each individual PCR product was sequenced using ABI BigDye Terminator Version 3.1 cycle sequencing kit (Applied Biosystems, Foster City, CA) on an ABI 3730 at the Center for Reproductive Biology Molecular Biology and Genomics Core facility at Washington State University. The obtained sequences were analyzed with BioEdit 7.0.1 software (Ibis Biosciences, Carlsbad, CA). The occurrence of mutations was evaluated according to the inspection of sequencing chromatographs, as containing one or both alleles. Each sample was sequenced three times with independently prepared genomic DNAs.
RNA extraction, cDNA synthesis and qRT-PCR. Total RNA from 100 spider mites per population was extracted using TRIZOL reagent (Invitrogen) following manufacturer's protocol. The quality of total RNA was checked by gel electrophoresis and spectrometry analyses. The total RNA was treated with DNase I (Ambion Inc., Austin, TX) to remove contaminating DNA. DNase I treated total RNA was used as a template for cDNA synthesizes by M-MLV reverse transcriptase (Promega, Madison, WI). qRT-PCR was performed using a CFX96 ™ Real-Time PCR Detection System (Bio-Rad Laboratories, Hercules, CA). Each qRT-PCR reaction (10 μl final volume) contained 5 μl iQ ™ SYBR Green Supermix (Bio-Rad Laboratories, Hercules, CA), 1.0 μl of cDNA, 3.6 μl ddH 2 O, and 0.4 μl forward and reverse gene specific primers (Table S4, stock 10 μM). An initial incubation of 95 °C for 3 min, followed by 40 cycles of 95 °C for 10 s, 55 °C for 60 s settings were used. The qRT-PCR for each sample was conducted with two technique replicates and three biological replicates. The no-template control and internal controls were included in each plate. Actin and rp49 were used as reference genes for internal controls 49 . Relative expression levels for target genes, in relation to two reference genes, actin and rp49 were calculated by the 2 −∆∆CT method 69 . Both the PCR efficiency and R 2 (correlation coefficient) value were taken into consideration in estimating relative quantities. PCR efficiency between 95% and 105% and R 2 value > 0.99 for each gene were considered as qualified for further analysis.