Skip to main content
  • Research article
  • Open access
  • Published:

Transcriptomic analysis of s-methoprene resistance in the lesser grain borer, Rhyzopertha dominica, and evaluation of piperonyl butoxide as a resistance breaker



The lesser grain borer, Rhyzopertha dominica is a serious pest of stored grains. Fumigation and contact insecticides play a major role in managing this pest globally. While insects are developing genetic resistance to chemicals, hormonal analogues such as s-methoprene play a key role in reducing general pest pressure as well as managing pest populations that are resistant to fumigants and neurotoxic contact insecticides. However, resistance to s-methoprene has been reported in R. dominica with some reports showing a remarkable high resistance, questioning the use of this compound and other related analogues in grain protection. The current study attempts to identify possible molecular mechanisms that contribute in resistance to s-methoprene in R. dominica.


Transcriptome analysis of resistant and susceptible strains of this pest species identified a set of differentially expressed genes related to cytochrome P450s, indicating their potential role in resistance to s-methoprene. Laboratory bioassays were performed with s-methoprene treated wheat grains in presence and absence of piperonyl butoxide (PBO), a cytochrome P450 inhibitor. The results indicate that PBO, when applied alone, at least at the concentration tested here, had no effect on R. dominica adult emergence, but has a clear synergistic effect to s-methoprene. The number of produced progeny decreased in presence of the inhibitor, especially in the resistant strain. In addition, we also identified CYP complement (CYPome) of R. dominica, annotated and analysed phylogenetically, to understand the evolutionary relationships with other species.


The information generated in current study suggest that PBO can effectively be used to break resistance to s-methoprene in R. dominica.


The lesser grain borer, Rhyzopertha dominica (F.) (Coleoptera: Bostrychidae) is among the most destructive pests of stored grains, with global distribution [1]. It is a primary feeder and infests a variety of stored products and related commodities [2], which are essential for human nutrition and global food security [1, 3]. Moreover, it is a primary colonizer, thus larvae and adults can easily penetrate the kernels even at low moisture content and complete their life cycle in intact whole grain kernels [2,3,4]. As a result, most life stages, especially the larvae, are unaffected by contact insecticides that are applied on the external part of the grain kernel [1]. Crucially, R. dominica has a rapid population growth resulting in devastating infestation levels, especially at optimal temperatures [1, 5]. Management of R. dominica in stored grain and other commodities have been investigated around the globe [1, 6]. In general, its control is currently based on two broad categories of insecticides, the fumigants [7] and contact insecticides [8]. However, it is now well-established that strains of R. dominica have developed resistance to both chemical and non-chemical treatments. In particular, high levels of resistance to phosphine [9,10,11], pirimiphos-methyl [12] and deltamethrin [7, 13] have been reported in many parts of the world, such as Australia, USA and Brazil [9,10,11]. At the same time, this species cannot be easily controlled by some “traditional” contact insecticides that are applied directly on grains, such as the organophosphorous compound pirimiphos-methyl [12] and the pyrethroid deltamethrin [7, 13]. Moreover, it is well-established that R. dominica is less susceptible than other major stored product insect species to non-chemical control methods, such as diatomaceous earths [14], which poses serious challenges to grain industry towards management of this species. Therefore, there is a demand to identify newer, reduced risk compounds that can be effectively used in controlling this important pest.

One of the newer active ingredients that have been registered in many countries for the control of R. dominica is the juvenile hormone analogue (JHA), s-methoprene, [15]. JHAs target and disrupt the endocrine system of insects by causing abnormal larval-pupal or nymphal-pupal development and/or even death [16]. In general, s-methoprene has many desirable characteristics, such as good environmental profile and extremely low mammalian toxicity [17, 18] and it is currently considered as a good alternative to many other conventional contact insecticides [15, 19,20,21,22]. It also exhibits a considerable residual efficacy on stored grains, thus holding a high potential as a grain protectant for long-term treatment [15, 23].

Although resistance to JHAs is not that frequent, resistance to pyriproxifen in the house fly Musca domestica L. (Diptera: Muscidae) and the whitefly Bemisia tabaci (Gennadius) (Hemiptera: Aleyrodidae) [24], as well as s-methoprene in mosquitoes [16] have been reported, suggesting that resistance may develop in the case of other species, including R. dominica. An s-methoprene resistant strain of R. dominica required a very high dose (40 mg kg− 1) for its control in wheat grain [25]. This dose rate is approximately 67 times higher than the registered rate applied in Australia, questioning the usage of this insecticide as a grain protectant. Moreover, resistance to s-methoprene may jeopardize the resistance management strategies to phosphine and neurotoxic insecticides [26], on which the inclusion of a JHA, e.g. on a rotation basis, is a key element.

Piperonyl butoxide (PBO), has been used extensively either alone or in combination with other active ingredients as a synergist in crop protection, especially to break resistance to specific group of insecticides such as pyrethroids that exhibits toxicity through mixed function oxidases including CYPs [27]. Several studies reported the interaction of PBO with cytochrome P450s [27, 28]. In the case of stored product protection, PBO has been successfully applied in many different cases [29,30,31].

The molecular mechanism of s-methoprene resistance has not been fully elucidated yet. In the vinegar fly, Drosophila melanogaster Meigen (Diptera: Drosophilidae), the absence of a so-called methoprene tolerant (MET) gene results in s-methoprene resistance [32, 33]. The protein (MET) encoded by the MET gene belongs to the family of basic helix-loop-helix (bHLH)-PAS transcriptional regulators that bind JH with high affinity [34]. MET forms homodimers (Gce in D. melanogaster forming heterodimer) in absence of ligand, i.e. Juvenile hormone III (JH-III), the growth juvenile hormone synthesized in most insects, or a synthetic mimic. In presence of either ligand, MET homodimer dissociates and their presence leads to dissociation of the MET dimer and thus binding with the ligand (JH-III or synthetic mimic). Ligand binding and immunoprecipitation assays where both MET monomers carry the V297F mutation, indicated resistance to s-methoprene thus they were not dissociated compared to the wild type counterpart [34]. Further experiments indicated that methoprene binds to PAS-B domain of the MET protein. Also, functional assays by knocking down MET in T. castaneum, render the insects resistant to the natural JH and as well as s-methoprene [35]. Alternatively, resistance to s-methoprene in other species has been associated with high activity of P450 monooxygenases and esterases, which probably also contribute to resistance to s-methoprene and other JHAs [36, 37]. However, detailed research revealing the exact relationship between s-methoprene and CYPs is not established, but it has been shown that P450s can metabolize JHAs, as in the case of pyriproxifen [38], which consists an indication that the same phenomenon may occur in the case of s-methoprene.

Resistance to s-methoprene has not been analysed yet in R. dominica, largely due to the lack of genomic resources for this pest species. RNA sequencing technologies have evolved rapidly in the last years [39]. They allow the study of transcriptomes without necessarily relying on a reference genome, thus greatly facilitating the study of several non-model species. Subsequently, comparison of gene transcription levels between insecticide resistant and insecticide-susceptible insect strains can lead to candidate genes that could play a role in the observed resistant phenotype. Such analysis has been performed in several insects and mites [40,41,42,43], providing not only a better understanding of insecticide resistance, but also valuable genomic resources that prove useful for studying different aspects of the biology of arthropods that constitute the most diverse animal clade [44,45,46].

In this regard, the aim of the present work was to investigate, for the first time, the mechanisms underlying s-methoprene resistance in R. dominica. We used s-methoprene-resistant and susceptible strains and compared their response to s-methoprene alone, but also in combination with PBO and mortality and progeny production were measured. The bioassays showed that the combined use of s-methoprene + PBO increased the efficacy of the former, thereby suggesting a possible involvement of CYPs in the resistance mechanism. Subsequently, we sequenced the transcriptomes of s-methoprene-resistant and susceptible strains and identified the Cytochrome P450 (CYP) genes. Interestingly, their analysis revealed that a number of them were significantly up-regulated in the resistant strain and are thus worth of further investigation to determine their role in insecticide resistance to JHAs.


Laboratory bioassays

Treatment effects were significant (Table 1). Parental mortality was low for 7, 14 and 21 days for both strains. Parental mortality for the control Lab-S was 0.1 and 12% for the Met-R. Moreover, for Lab-S and Met-R the lowest parental mortality was 6.7 and 2.2 and the highest 26.7 and 17.8 respectively (Additional file 1: Fig. S1). Regarding progeny production counts, adult emergence was generally higher in the case of the resistant strain, as compared to the respective figures of the susceptible strain, even in the untreated grains (Fig. 1). Moreover, the application of PBO alone, for both strains, had no effect, as the numbers of adults that had been emerged after the termination of the incubation period were extremely high (> 150 adults/vial), and comparable to those in the controls (Fig. 1). Still, for the resistant strain, the application of PBO alone caused a slight reduction in progeny production, in comparison with the control vials. In the case of the Lab-S strain, the combination of PBO with s-methoprene gave similar results with the application of s-methoprene alone. For this strain, when s-methoprene was applied either alone or in combination, progeny production was generally higher at 0.01 mg/kg than that for the other concentrations. Nevertheless, for the susceptible strain, progeny production ranged between 0.2 and 2.3 adults/vial (Fig. 1). In contrast, for the resistant strain, when s-methoprene was applied alone, progeny production was significantly lower than that in the control vials (Fig. 1). However, there was a considerably high offspring emergence, regardless of the concentration. The increase of the concentration from 1 to 30 mg/kg resulted in a gradual decrease on the number of emerged R. dominica adults, from 122 to 33 individuals/vial. Similarly, when s-methoprene was applied with PBO, the increase in the concentrations reduced progeny production from 120 to 19 individuals/vial (Fig. 1). Furthermore, for the two lowest s-methoprene concentrations, progeny production was not affected, regardless of the presence of PBO. Nevertheless, for the two higher concentrations, progeny production of R. dominica was considerably lower when s-methoprene was applied in combination with PBO, than for the application of s-methoprene alone (Fig. 1).

Table 1 ANOVA parameters for progeny production of R. dominica susceptible (Lab-S) and resistant strain (Met-R) (error df=80)
Fig. 1
figure 1

Mean number (±SEM) of Rhyzopertha dominica progeny production (expressed as adults/vial) for the susceptible (Lab-S) and resistant (Met-R) strain for all the combinations tested (control, 0.01 mg/kg, 0.03 mg/kg, 0.1 mg/kg, 0.3 mg/kg, PBO, 0.01 mg/kg + PBO, 0.03 mg/kg + PBO, 0.1 mg/kg + PBO, 0.3 mg/kg + PBO for susceptible and control, 1 mg/kg, 3 mg/kg, 10 mg/kg, 30 mg/kg, PBO, 1 mg/kg + PBO, 3 mg/kg + PBO, 10 mg/kg + PBO, 30 mg/kg + PBO for resistant). Within each bar and strain, means followed by the same lowercase letter do not differ significantly according to Tukey Kramer HSD test at P< 0.05. Where no letter exist, no significant differences were noted. Means with asterisk (*) for the application with s-methoprene alone are significantly different for the respective mean of the combination with s-methoprene and PBO at the resistant strain (Met-R)

Transcriptome sequencing

In order to better study the molecular basis of the observed resistance, the transcriptome of R. dominica was sequenced, yielding a total of > 688 million Illumina reads. These reads were then assembled de novo with Trinity since there is no available reference genome sequence. The assembled transcriptome contained a total of 117,265 putative transcripts (Table 2). The quality of the assembly is very good, as evidenced by the BUSCO analysis [47], which showed that > 98% of the conserved insect genes are present in the assembly (Table 3).

Table 2 Transcriptome assembly summary
Table 3 Detailed RNA sequencing results for each R. dominica strain

After calculating transcript abundance a Principal Component Analysis (PCA) was run in order to verify the quality of the biological replicates. It is evident that the replicates of the same strain clustered together, but also are separated from the replicates of the other strain (Additional file 2: Fig. S2). These results show that the sequencing data are of good quality and can be used in downstream analyses.

Investigating target site-mediated resistance

The sequence polymorphism analysis as well as expression levels in the MET gene between the Lab-S and Met-R R. dominica strains did not detect any significant differential expression. However, examination of the open reading frame (ORF) of MET between the two strains revealed the occurrence of a non-synonymous amino acid substitution at position 489 of the amino-acid sequence in the Met-R strain. The observed substitution leads to the replacement of a Pro by Leu. However, this mutation is not fixed in Met-R, it is present in only 33% of the reads, and, finally, is located outside of the PAS-B conserved domain.

Investigation of non-target site resistance mechanisms based on differential expression and qPCR validation

Differential expression (DE) analysis was done on all the 117,265 assembled transcripts, at the unigene level. This analysis showed that 275 unigenes were up-regulated in the Met-R strain compared to Lab-S, whereas another 190 were down-regulated (Fig. 2). No significantly over-represented GO terms or KEGG pathways were found in either the up-or down-regulated set of genes (padj < 0.01).

Fig. 2
figure 2

Overview of the differentially expressed (|log2FC| > 2 and also p-value < 0.001) genes between the resistant and the susceptible to s-methoprene strains of R. dominica. In total, there were 465 differentially expressed unigenes, of which 276 are up-regulated in the resistant strain, whereas the remaining 190 are up-regulated in the susceptible strain. The data points corresponding to P450s have been colored as red, whereas the one corresponding to the UGT is colored as purple

Interestingly, we identified a number of up or down-regulated unigenes that have a similarity to detoxification enzymes (Table 5). These include six CYPs (DN26728_c0_g1, DN29475_c1_g7, DN28703_c3_g1, DN23343_c0_g1, DN28703_c3_g3, DN26679_c1_g1), one glutathione S-transferase (GST) (TRINITY_DN20738_c0_g1), and one UDP-glucosyltransferase (UGT) (DN28972_c1_g2). The CYPs as well as the UGT were up-regulated in the Met-R strain, whereas the GST was up-regulated in the Lab-S strain. The difference in expression levels was statistically significant for all these unigenes (FDR < 0.05). The over-expression levels of the identified CYPs were validated by qPCR with CYP6BQ11 (DN26728_c0_g1), CYP6RU (DN28703_c3_g1 and DN23 343_c0_g1) and CYP3747A (DN26679_c1_g1) displaying significant (p=value < 0.05) up-regulation of > 10-, 4- and 3-fold in the Met-R strain, compared to the Lab-S strain (Fig. 3).

Fig. 3
figure 3

Relative expression levels of the six CYPs. Expression levels are depicted relative to Lab-S reference susceptible strains. Error bars represent 95% confidence intervals. Asterisks indicate significantly different expression (p-value < 0.05)

Detailed study of putative CYPs

Rhyzopertha dominica transcripts containing the InterPro domain IPR001128, were searched and annotated as putative CYPs or CYP fragments. The analysis revealed 396 probable CYP isoforms of R. dominica putatively originating from 111 unigenes. Maximum Likelihood phylogenetic analysis was performed on the largest isoform from each unigene, using the T. castaneum CYP genes [48] as a reference. All the R. dominica CYPs were classified into one of the four known CYP clans existing in T. castaneum (Fig. 4, Table 4). Furthermore, this analysis revealed at least eight R. dominica-specific clades in Clans 3 and 4 for some of which a clear classification within the respective clade was not possible. In addition, the phylogenetic analysis also shows that there are four different unigenes in R. dominica that cluster with the T. castaneum CYP12H1. This is an indication of probable duplication events that led to multiple copies of this CYP gene in R. dominica.

Fig. 4
figure 4

Phylogenetic analysis of the CYP genes identified in R. dominica. This analysis showed that all identified R. dominica CYPs could be classified into one of the known T. castaneum clans. Furthermore, the differentially expressed CYPs belong to Clan 3 (four unigenes) and Clan 4 (two unigenes). All R. dominica genes were classified into one of the four known CYP clans previously found in the beetle T. castaneum. Bootstrap values > 75% are represented as black dots, whereas nodes with bootstrap support between 50 and 75% are shown as grey dots. Nodes with bootstrap support < 50% are collapsed. The R. dominica-specific expansions in Clans 3 and 4 containing the up-regulated CYPs are highlighted in light orange and light green, respectively. CYPs whose log2FC is > 2 are marked with a red asterisk, whereas those with a log2FC between 1 and 2 are marked with a red triangle. The scale bar is in substitutions per site

Table 4 Summary of the phylogenetic analysis of R. dominica P450s

Interestingly, two of the identified CYPs were significantly up-regulated (FDR < 0.001, log2FC > 2) in the Met-R strain. Another four also appear to be significantly up-regulated, albeit at a lower degree (FDR < 0.05, log2FC > 1.44). Four of the six CYPs belong to Clan 3, whereas the other two belong to Clan 4 (Fig. 3, Table 5). A more precise placement of these clades was not possible due to the low bootstrap support values (< 50%) of the respective branches. Nevertheless, expert manual curation by Dr. David Nelson annotated these genes as similar to the genes CYP6BQ11, CYP3747A (from Oryctes borboronicus), CYP6RU (from Photinus pyralis), and CYP6RU1 (from Photinus pyralis) (Table 5; Additional file 3: Table S1).

Table 5 Up-regulated detoxification enzymes in the Met-R strain


Τhe frequency of cases of insecticide resistance of insects infesting stored products has been increased in the last decades [12, 26, 49,50,51]. S-methoprene is an insect growth regulator which plays a pivotal role in mitigating resistance to several contact insecticides and fumigants [15, 52]. Although it has a unique mode of action, and it has not been previously used as grain protectant, there are reports of high levels of resistance to s-methoprene in R. dominica [24, 53], which may question its use in the near future [26]. While many studies have focused on the phenotypic characterization of resistance, the current work elucidates molecular mechanisms of resistance to s-methoprene in R. dominica, with a perspective of developing suitable resistance management practices.

Our study clearly indicated that the simultaneous application of s-methoprene + PBO, increased the insecticidal effect of s-methoprene (Fig. 1) and all the above clearly indicate that PBO, which when applied alone, at least at the concentration tested here, had no effect on R. dominica adult emergence, but has a clear synergistic effect to s-methoprene.

The use of PBO as a synergist has been extensively used in stored product protection, but most of the studies available are about pyrethroids. For instance, the application of PBO with natural pyrethrum were found to increase the efficacy of diatomaceous earths for the control of R. dominica on different grains [31]. Similar results have been reported for the application of natural pyrethrum alone [54]. Deltamethrin resistance has been shown to reduce from 223-fold to 21-fold using the CYP inhibitor PBO against a pyrethroid resistant population of the cotton armyworm, Helicoverpa armigera (Hübner) (Lepidoptera: Noctuidae) [55]. In the case of stored product insects, the granary weevil, Sitophilus granarius (L.) (Coleoptera: Curculionidae) was tested with PBO and fenitrothion and it was found that there is a positive synergism between them [56]. Our results suggest that this combination can also be effective in the case of resistance to JHA by stored product insects, but, to our knowledge such an approach has not been implemented yet.

Sequence analysis of the MET gene identified a P489L substitution in the resistant Met-R strain, but not in the susceptible Lab-S. A mutational change at position 297 in the MET protein was reported earlier in s-methoprene-resistant T. castaneum that has explicitly exhibited reduced binding affinity to s-methoprene [34]. However, the herein identified P489L mutation is located at the C-terminus of the gene and outside of the PAS-B domain that has been previously implicated in ligand binding. The functional role of P489L and its contribution to resistance remains to be investigated.

Our analysis identified a total of 111 R. dominica CYPs, a number relatively close to the 143 CYPs identified in the closely related T. castaneum [48]. The missing genes in R. dominica are most probably due to the fact that not all the CYP genes were transcribed at the samples we sequenced. Six of these CYPs were up-regulated in the resistant Met-R strain (Fig. 2), with four of them belonging to Clan 3 and the other two to Clan 4. However, their comparison to T. castaneum CYPs showed that they cannot be reliably grouped to any family within these clans (Fig. 3). Evidence supports that CYPs metabolize JHAs, such as in the case of pyriproxifen [38]. Additionally, heterologously expressed CYPs from A. gambiae were shown to be capable of metabolizing pyriproxifen, with CYP6P3 to be the strongest metabolizer [57]. Moreover, microsomal CYPs are capable of metabolizing s-methoprene when incubated with housefly microsomes [58]. Application of sub-lethal concentrations of s-methoprene in Sf9 cells and the fall armyworm indicated induction of the expression of CYPs in Sf9 cells, most of which belong to the CYP9 family, whereas in live insects CYP9A28 was differentially expressed in response to s-methoprene [59]. Thus, CYPs could be potentially involved in s-methoprene resistance and are therefore worth of a more detailed analysis. To this end, and in order to be able to properly classify each of the identified R. dominica CYPs, we conducted a Maximum Likelihood phylogenetic analysis with the manually curated set of CYPs of T. castaneum [48]. The results of this analysis showed that all six up-regulated CYPs belong to three well-supported clades which, nevertheless, only contained R. dominica genes (Fig. 3). This, of course, did not allow for a more precise placement of these clades within the respective CYP clans. Nevertheless, expert manual annotation by Dr. David Nelson assigned specific functions to each one of them, thus providing hints for their possible function (Additional file 4: Table S2).

More specifically, CYPs belonging to the CYP6 family (Clan 3) have been characterized and shown to metabolize xenobiotics and plant allelochemicals in several insect species, such as M. domestica, the African malaria mosquito, Anopheles gambiae Giles (Diptera: Culicidae) and others [60]. However, there was no detection of a specific CYP gene capable of metabolizing s-methoprene. For example, the M. domestica CYP6A1 was not able to metabolize s-methoprene and hydroprene, while other JHs such as juvenile hormones I and III were metabolized [61]. Nevertheless, transcriptomic analysis of a pyriproxyfen (another JHA insecticide)-resistant strain of the greenhouse whitefly, Trialeurodes vaporariorum Westwood (Hemiptera: Aleurodidae) suggested that the most highly up-regulated CYPs (log2FC between 2.68 and 2.91) also belonged to the CYP6 family [62], but qPCR analysis indicated that a CYP belonging to the CYP4 family (Clan 4) is highly upregulated in the pyriproxyfen selected strain. Here, the expression levels of 6 CYPs from R. dominica were validated by qPCR indicating two of them to be significantly over-expressed in the Met-R strain in comparison to the Lab-S. Given the synergistic effect of PBO that argues in favor of a P450-mediated resistance and validated over-expression of four CYPs, point towards their functional expression and characterization which will give evidence for their role in metabolic activity and implication in the observed resistance.


The results of the present study indicate that resistance to s-methoprene is potentially mediated by cytochrome CYPs. Moreover, our bioassay data suggested that the simultaneous application of PBO and s-methoprene can be used with success to mitigate resistance to s-methoprene in R. dominica. In order to investigate whether CYPs are involved in resistance we sequenced a susceptible (Lab-S) and a resistant (Met-R) to s-methoprene strain of R. dominica, using RNA-seq. Data analysis identified a number of up-regulated genes that bear significant similarity with CYPs and could thus be involved in the detoxification of s-methoprene. Most importantly, the herein generated transcriptome assembly is the only genomic resource for R. dominica and it can serve as a reference in future projects studying the biology of this important pest. Additionally, our results suggest that PBO acts as a “resistance breaker” and should therefore be considered towards the direction of resistance management. This is particularly important in the case of s-methoprene, as it is classified among the insecticides with the lowest mammalian toxicity that are currently in use as grain protectants.


Reagents and chemicals

The pure analytical grade chemicals of s-methoprene and PBO-8 were obtained from Dow Agrosciences Ltd. (CPC2 Capital Park, Fulbourn, Cambridge, England, CB21 5XE). In addition, the commercial formulations of these two chemicals, such as Diacon IGRTR and PBO8-Synergist, containing 33.6% active ingredient (a.i) of s-methoprene and 91.3% of PBO, respectively, were used for the tests.

Insect strains

Two strains of R. dominica, which are susceptible (Lab-S) and resistant to s-methoprene (Met-R), respectively, were used in this study. The Lab-S was collected form a grain storage shed at Oakey in 1971, whereas, Met-R strain was collected from Roma, Queensland [26]. Since then, the insect populations of these two strains were reared on whole wheat kernels and maintained at standard room temperature and relative humidity (RH) at Queensland Department of Agriculture and Fisheries (QDAF), Australia. Adult beetles less than three weeks old were randomly selected from the culturing jars and used in the bioassays.

Laboratory bioassays

Two sets of s-methoprene concentrations were used in the bioassays. These include, 0, 0.01, 0.03, 0.1 and 0.3 mg kg− 1 for Lab-S and 0, 1, 3, 10 and 30 mg kg− 1 for Met-R. For PBO bioassay, the wheat grains were applied with a recommended label rate for combinations, 0.013 lt per 45.3 kg of wheat. Untreated clean and infestation-free wheat grains were used for treatments. The moisture content of the grain was adjusted to 13.5% before the initiation of the experiment. Three lots of wheat containing 2 kg each were sprayed with different dose rates of s-methoprene alone, PBO alone and the combination of s-methoprene + PBO; hence, the combinations of the formulations were: a) s-methoprene alone, b) PBO alone and c) s-methoprene with PBO, in all possible combinations (control, 0.01 mg/kg, 0.03 mg/kg, 0.1 mg/kg, 0.3 mg/kg, PBO, 0.01 mg/kg + PBO, 0.03 mg/kg + PBO, 0.1 mg/kg + PBO, 0.3 mg/kg + PBO for Lab-S and control, 1 mg/kg, 3 mg/kg, 10 mg/kg, 30 mg/kg, PBO, 1 mg/kg + PBO, 3 mg/kg + PBO, 10 mg/kg + PBO, 30 mg/kg + PBO for Met-R). The required volume of each treatment including the combinations, was applied using a specialized airbrush (Badger 100, Kyoto BD-183 K Grapho-tech, Japan). An additional Series of 2 kg wheat lots, sprayed with water in parallel to each treatment, was used as untreated control. Twenty grams (20 g) of grain samples from each treatment was selected for bioassays. This sample was placed inside cylindrical bioassay vials (3 cm in diameter, 8 cm in height). Ten adults of R. dominica were released into each vial. The vials were then placed in incubators set at 27.5 °C and 75% RH. The mortality of adults was recorded after 7, 14 and 21 days of exposure. Thereafter, all parental adults were removed from the vials, and the vials with the treated grain were returned to the same incubators and maintained for 65 d more to ensure that the immatures in the treated grain will develop up to the adult stage. Then, the number of adults that emerged in treatments and control were compared and per cent reduction in progeny production was estimated. The entire experiment was repeated three times (jars) with each containing 3 subreplicates (vials).

RNA isolation, library construction sequencing and qPCR validation

Ten 2nd to 3rd instar larvae of R. dominica of Lab-S and Met-R strains were pooled respectively and preserved in RNA later, and total RNAs of each was extracted using the GeneJet RNA Purification kit (ThermoScientific), according to the manufacturer’s protocol. The extracted RNA was treated with Turbo DNase (Ambion), in order to remove any traces of genomic DNA. The purity and concentration of RNA were estimated using a Nanodrop spectrophotometer based on 260/280 and 260/230. RNA samples were sent to Macrogen (Korea) for mRNA paired-end library construction with the Illumina Truseq stranded mRNA sample preparation kit, following the manufacturer’s instructions. Each library was sequenced with the paired-end method for a read length of 100 bp. Two μg of RNA was used for cDNA synthesis using the reverse transcriptase kit from Minotech (Heraklion, Greece), according to the manufacturer’s instructions. qPCR validation was conducted for a subset of genes. The primers used are shown on Additional file 6: Table S4. Briefly, a 5-fold dilution series of pooled cDNA was used to assess the efficiency of the qPCR reaction for each gene-specific primer pair. A no template control (NTC) was also included to detect possible contamination. The reactions consisted of 0.6 μM primers each, and Kapa SYBR FAST qPCR Master Mix (Kapa-Biosystems). Experiments were performed using 3 biological and 3 technical replicates for each gene. The levels of the validated genes were measured by Real-time qPCR (RT-qPCR) amplification on a CFXConnect (BioRad). Relative expression levels were calculated as previously described [63].

Computational analyses

RNAseq reads from both strains (total of ~ 688 million reads) were assembled with Trinity v2.5.1 [64], using parameters “--seqType fq --SS_lib_type RF --max_memory 350G --CPU 24”. InterProScan v5.28–67 [65] was used in order to identify conserved domains within each assembled transcript. Moreover, BLAST v2.8.0+ [66] searches were run in order to identify similarities using the Uniref50 database that is specifically built for similarity-based functional annotation [67].

Transcript abundance was estimated with Kallisto [68]. Next, the scripts bundled with Trinity were used for running the differential expression analysis with EdgeR [69] in order to find transcripts that were differentially expressed between the resistant and the susceptible strain (FDR < 0.05). Custom Perl and bash scripts were used for parsing the EdgeR output and identifying genes of interest. Gene Ontology (GO) term analyses were done using gProfiler [70].

For the detection of polymorphisms in the methoprene tolerant gene we firstly mapped the raw reads to the Trinity transcripts using hisat2 [71], then generated a mpileup file with samtools [72], and searched for SNPs with VarScan v2.4.4 [73]. Finally, the identified SNPs were visually inspected across the extracted methoprene tolerant transcript using samtools and the data were loaded into the Integrative Genomics Viewer v2.6.3 [74].

In order to identify transcripts with similarity to cytochrome P450 (CYP) genes we first ran the TransDecoder program that is bundled with Trinity v2.8.5 [64] and obtained the encoded peptides in each transcript. Subsequently, putative CYP-related proteins were identified by the presence of the IPR001128 InterPro domain, in the InterProScan output file. The curated CYPs set identified in T. castaneum were obtained from [48] and used as a reference for classifying the herein identified R. dominica CYPs. Finally, the early-diverged CYP51A1 [75] from Homo sapiens was used as an outgroup. Multiple sequence alignment (Additional file 5: Table S3) was performed with MAFFT v7.271 [76] with parameters “--auto --threads 8” and trimming was done with Trimal v1.2rev59 [77], with parameters “--gt 0.50”. A Maximum Likelihood phylogeny with 100 bootstrap replicates was inferred with RAxML v8.2.11 [78], with parameters “-m PROTGAMMAAUTO”. Branches with < 50% bootstrap support were collapsed with TreeGraph2 [79] and the resulting Newick tree was loaded to a locally deployed instance of EvolView v2 [80] for post-processing. The vector graphics editor Inkscape v0.92 was used for the final polishing.

Bioassay data analysis

The data of progeny production were analyzed separately for each strain using ANOVA to test the treatment effects. When preliminary tests indicated that variances were not equal, the data were transformed to log (x+ 1) (for the susceptible strain O’Brien test: F=1.01, P=0.437; for the resistant strain O’Brien test: F=1.84, P=0.073). Means were separated by using the Tukey-Kramer HSD test at the 5% level. For each strain the Student’s t-test was used to determine differences between s-methoprene alone and s-methoprene with PBO. Statistical analysis was performed by using the JMP 7 software (SAS Institute Inc., Cary, NC, USA).

Availability of data and materials

The sequencing reads are available from the Sequence Read Archive (SRA) under the bioproject accession PRJNA605183. Additional data and custom scripts are available upon request.



Piperonyl butoxide


Juvenile hormone analogue


Differential expression


  1. Edde PA. A review of the biology and control of Rhyzopertha dominica (F.) the lesser grain borer. J Stored Prod Res. 2012;48:1–18.

    Article  Google Scholar 

  2. Hagstrum DW, Subramanyam B. Stored-product insect resource. In: AACC International, saint 943 Paul. USA: MN; 2009.

    Google Scholar 

  3. Kavallieratos NG, Athanassiou CG, Arthur FH, Throne JE. Lesser grain borers, Rhyzopertha dominica, select rough rice kernels with cracked hulls for reproduction. J Insect Sci. 2012;12:38.

    Article  PubMed  PubMed Central  Google Scholar 

  4. Rahim M. Biological activity of azadirachtin-enriched neem kernel extracts against Rhyzopertha dominica (F.) (Coleoptera: Bostrichidae) in wheat. J Stored Prod Res. 1998;34:123–8.

    Article  CAS  Google Scholar 

  5. Sakka M, Athanassiou CG. Competition of three stored-product bostrychids on different temperatures and commodities. J Stored Prod Res. 2018;79:34–9.

    Article  Google Scholar 

  6. Athanassiou CG, Arthur FH. Bacterial Insecticides and Inert Materials. In: Athanassiou C, Arthur F, editors. Recent advances in stored product protection. Berlin, Heidelberg: Springer; 2018.

    Chapter  Google Scholar 

  7. Daglish GJ, Nayak MK. Prevalence of resistance to deltamethrin in Rhyzopertha dominica (F.) in eastern Australia. J Stored Prod Res. 2018;78:45–9.

    Article  Google Scholar 

  8. Chen CY, Chen MF. Susceptibility of field populations of the lesser grain borer, Rhyzopertha dominica (F.), to deltamethrin and spinosad on paddy rice in Taiwan. J Stored Prod Res. 2013;55:124–7.

    Article  Google Scholar 

  9. Collins PJ, Daglish GJ, Bengston M, Lambkin TM, Pavic H. Genetics of resistance to phosphine in Rhyzopertha dominica (Coleoptera: Bostrichidae). J Econ Entomol. 2002;95:862–9.

    Article  PubMed  Google Scholar 

  10. Lorini I, Collins PJ, Daglish GJ, Nayak MK, Pavic H. Detection and characterisation of strong resistance to phosphine in Brazilian Rhyzopertha dominica (F.) (Coleoptera: Bostrychidae). Pest Manag Sci. 2007;63:358–64.

    Article  CAS  PubMed  Google Scholar 

  11. Opit GP, Phillips TW, Aikins MJ, Hasan MM. Phosphine resistance in Tribolium castaneum and Rhyzopertha dominica from stored wheat in Oklahoma. J Econ Entomol. 2012;105:1107–14.

    Article  CAS  PubMed  Google Scholar 

  12. Guedes PNC, Dover BA, Kambhamoatik S. Resistance to chlorpyrifos-methyl, pirimiphos-methyl, and malathion in Brazilian and U.S. populations of Rhyzopertha dominica (Coleoptera: Bostrichidae). J Econ Entomol. 1996;89:27–32.

    Article  CAS  Google Scholar 

  13. Lorini I, Galley DJ. Deltamethrin resistance in Rhyzopertha dominica (F.) (Coleoptera: Bostrichidae), a pest of stored grain in Brazil. J Stored Prod Res. 1999;35:37–45.

    Article  CAS  Google Scholar 

  14. Fields P, Korunic Z. The effect of grain moisture content and temperature on the efficacy of diatomaceous earths from different geographical locations against stored-product beetles. J Stored Prod Res. 2000;36:1–13.

    Article  Google Scholar 

  15. Wijayaratne LKW, Arthur FH, Whyard S. Methoprene and control of stored-product insects. J Stored Prod Res. 2018;76:161–9.

    Article  Google Scholar 

  16. Kamita SG, Samra AI, Liu JY, Cornel AJ, Hammock BD. Juvenile hormone (JH) esterase of the mosquito Culex quinquefasciatus is not a target of the JH analog insecticide methoprene. PLoS One. 2011;6:e28392.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  17. Hawkins DR, Weston KT, Chasseaud LF, Franklin ER. Fate of Methoprene (isopropyl (2E, 4E)-11-methoxy-3, 7, 11-trimethyl-2, 4-dodecadienoate) in rats. J Agric Food Chem. 1977;25:398–403.

    Article  CAS  PubMed  Google Scholar 

  18. Oberlander H, Silhacek DL, Shaaya E, Ishaaya I. Current status and future perspectives of the use of insect growth regulators for the control of stored product insects. J Stored Prod Res. 1997;33:1–6.

    Article  CAS  Google Scholar 

  19. Arthur FH. Evaluation of methoprene alone and in combination with diatomaceous earth to control Rhyzopertha dominica (Coleoptera: Bostrichidae) on stored wheat. J Stored Prod Res. 2004;40:485–98.

    Article  CAS  Google Scholar 

  20. Chanbang Y, Arthur FH, Wilde GE, Throne JE. Efficacy of diatomaceous earth and methoprene, alone and in combination, against Rhyzopertha dominica (F.) (Coleoptera: Bostrychidae). J Stored Prod Res. 2007;43:396–401.

    Article  CAS  Google Scholar 

  21. Chanbang Y, Arthur FH, Wilde GE, Throne JE, Subramanyam B. Susceptibility of eggs and adult fecundity of the lesser grain borer, Rhyzopertha dominica, exposed to methoprene. J Insect Sci. 2008;8:48–54.

    Article  PubMed Central  Google Scholar 

  22. Athanassiou CG, Arthur FH, Throne JE. Efficacy of layer treatment with methoprene for control of Rhyzopertha dominica (Coleoptera: Bostrychidae) on wheat, rice and maize. Pest Manag Sci. 2011;67:380–4.

    Article  CAS  PubMed  Google Scholar 

  23. Daglish GJ, Eelkema M, Harrison LM. Chlorpyrifos-methyl plus either methoprene or synergized phenothrin for control of Coleoptera in maize in Queensland, Australia. J Stored Prod Res. 1995;31:235–41.

    Article  CAS  Google Scholar 

  24. Palli SR. Recent advances in the Mode of Action of Juvenile Hormones and Their Analogs. In: Ishaaya I, Horowitz AR, editors. Biorational of Arthropod Pests: Application and Resistance Management; 2009. p. 111.

    Chapter  Google Scholar 

  25. Daglish GJ, Holloway JC, Nayak MK. Implications of methoprene resistance for managing Rhyzopertha dominica (F.) in stored grain. J Stored Prod Res. 2013;54:8–12.

    Article  CAS  Google Scholar 

  26. Collins PJ, Schlipalius DI. Insecticide Resistance. Athanassiou C, Arthur F (eds) Recent Advances in Stored Product Protection. Heidelberg: Springer, Berlin; 2018.

    Google Scholar 

  27. Hodgson E and Levi PE, 1998. Interactions of piperonyl butoxide with cytochrome P450, in Piperonyl Butoxide. Jones DJ (ed). The insecticide synergist. Academic Press, London, 1998; 41–53.

  28. Casida JE. Mixed-function oxidase involvement in the biochemistry of insecticide synergists. J Agric Food Chem. 1970;18:753–72.

    Article  CAS  PubMed  Google Scholar 

  29. Athanassiou CG, Kavallieratos NG, Andris NS. Insecti-cidal effect of three diatomaceous earth formulations against adults of Sitophilus oryzae (Coleoptera: Curculionidae) and Tribolium confusum (Coleoptera: Tenebrionidae) on oat, rye and triticale. J Econ Entomol. 2004;97:2160–7.

    Article  PubMed  Google Scholar 

  30. Kavallieratos NG, Athanassiou CG, Pashalidou FG, Andris NS, Tomanović Ž. Influence of grain type on the insecticidal efficacy of two diatomaceous earth formulations against Rhyzopertha dominica (F) (Coleoptera: Bostrychidae). Pest Manag Sci. 2005;61:660–6.

    Article  CAS  PubMed  Google Scholar 

  31. Athanassiou CG, Kavallieratos NG. Insecticidal effect and adherence of PyriSec in different grain commodities. Crop Prot. 2005;24:703–10.

    Article  CAS  Google Scholar 

  32. Wilson TG, Fabian J. A Drosophila melanogaster mutant resistant to a chemical analogue of juvenile hormone. Dev Biol. 1986;118:190–201.

    Article  CAS  PubMed  Google Scholar 

  33. Wilson TG, Ashok M. Insecticide resistance resulting from absence of target-site gene product. Proc Natl Acad Sci. 1998;95:14040–4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Charles JP, Iwema T, Epa VC, Takaki K, Rynes J, Jindra M. Ligand-binding properties of a juvenile hormone receptor. Methoprene-tolerant PNAS. 2011;108:21128–33.

    Article  CAS  PubMed  Google Scholar 

  35. Konopova B, Jindra M. Juvenile hormone resistance gene Methoprene-tolerant controls entry into metamorphosis in the beetle Tribolium castaneum. PNAS. 2007;104:10488–93.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Cerf DC, Georghiou GP. Cross resistance to juvenile hormone analogues in insecticide-resistant strains of Musca domestica. Pestic Sci. 1974;5:759–67.

    Article  CAS  Google Scholar 

  37. Brown TM, Brown AWA. Experimental in-duction of resistance to a juvenile hormone mimic. J Econ Entomol. 1974;67:799–801.

    Article  CAS  PubMed  Google Scholar 

  38. Zhang L, Kasai S, Shono T. In vitro metabolism of pyriproxifen by microsomes from susceptible and resistant housefly larvae. Arch Insect Biochem Physiol. 1998;37:215–24.

    Article  CAS  PubMed  Google Scholar 

  39. Stuart T, Satija R. Integrative single-cell analysis. Nat Rev Genet. 2019;20:257–72.

    Article  CAS  PubMed  Google Scholar 

  40. Crawford JE, Guelbeogo WM, Sanou A, Traoré A, Vernick KD, Sagnon NF, Lazzaro BP. Novo transcriptome sequencing in Anopheles funestus using illumine RNA-Seq technology. PLoS One. 2010;5:e14202.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  41. Karatolos N, Pauchet Y, Wilkinson P, Chauhan R, Denholm I, Gorman K, Nelson DR, Bass C, RHF C, Williamson MS. Pyrosequencing the transcriptome of the greenhouse whitefly, Trialeurodes vaporariorum reveals multiple transcripts encoding insecticide targets and detoxifying enzymes. BMC Genomics. 2011;12:56.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Bajda S, Dermauw W, Greenlagh R, Nauen R, Tirry L, Clark RM, Van Leeuwen T. Transcriptome profilimg of a spirodiclofen susceptible and resistant strain of the European red mite Panonychus ulmi using strand-specific RNA-seq. BMC Genomics. 2015;16:974.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  43. Ilias A, Lagnel J, Kapantaidaki DE, Roditakis E, Tsigenopoulos CS, Vontas J, Tsagkarakou A. Transcription analysis of neonicotinoid resistance in Mediterranean (MED) populations of B. tabaci reveal novel cytochrome P450s, but no nAChR mutations associated with the phenotype. BMC Genomics. 2015;16:1–23.

    Article  CAS  Google Scholar 

  44. i5K Consortium. The i5K Initiative: advancing arthropod genomics for knowledge, human health, agriculture, and the environment. J Hered. 2013;104:595–600.

    Article  PubMed Central  Google Scholar 

  45. Misof B, Liu S, Meusemann K, Peters RS, Donath A, Mayer C, Frandsen PB, Ware J, Flouri T, Beutel RG, Niehuis O, Petersen M, Izquierdo-Carrasco F, Wappler T, Rust J, Aberer AJ, Aspöck U, Aspöck H, Bartel D, Blanke A, Berger S, Böhm A, Buckley TR, Calcott B, Chen J, Friedrich F, Fukui M, Fujita M, Greve C, Grobe P, Gu S, Huang Y, Jermiin LS, Kawahara AY, Krogmann L, Kubiak M, Lanfear R, Letsch H, Li Y, Li Z, Li J, Lu H, Machida R, Mashimo Y, Kapli P, DD MK, Meng G, Nakagaki Y, Navarrete-Heredia JL, Ott M, Ou Y, Pass G, Podsiadlowski L, Pohl H, von Reumont BM, Schütte K, Sekiya K, Shimizu S, Slipinski A, Stamatakis A, Song W, Su X, Szucsich NU, Tan M, Tan X, Tang M, Tang J, Timelthaler G, Tomizuka S, Trautwein M, Tong X, Uchifune T, Walzl MG, Wiegmann BM, Wilbrandt J, Wipfler B, Wong TK, Wu Q, Wu G, Xie Y, Yang S, Yang Q, Yeates DK, Yoshizawa K, Zhang Q, Zhang R, Zhang W, Zhang Y, Zhao J, Zhou C, Zhou L, Ziesmann T, Zou S, Li Y, Xu X, Zhang Y, Yang H, Wang J, Wang J, Kjer KM, Zhou X. Phylogenomics resolves the timing and pattern of insect evolution. Science. 2014;346:763–7.

    Article  CAS  PubMed  Google Scholar 

  46. Dowling D, Pauli T, Donath A, Meusemann K, Podsiadlowski L, Petersen M, Peters RS, Mayer C, Liu S, Zhou X, Misof B, Niehuis O. Phylogenetic origin and diversification of RNAi pathway genes in insects. Genome Biol Evol. 2016;8:3784–93.

    CAS  PubMed  PubMed Central  Google Scholar 

  47. Waterhouse RM, Seppey M, Simão FA, Manni M, Ioannidis P, Klioutchnikov G, Kriventseva EV, Zdobnov EM. BUSCO applications from quality assessments to gene prediction and phylogenomics. Mol Biol Evol. 2017;35:543–8.

    Article  PubMed Central  CAS  Google Scholar 

  48. Zhu F, Moural TW, Shah K, Palli SR. Integrated analysis of cytochrome P450 gene superfamily in the red flour beetle, Tribolium castaneum. BMC Genomics. 2013;14:174.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  49. Champ BR, Campbell-Brown MJ. Insecticide resistance in Australian Tribolium castaneum (Herbst) (Coleoptera, Tenebrionidae) — II: Malathion resistance in eastern Australia. J Stored Prod Res. 1970;6:111–31.

    Article  CAS  Google Scholar 

  50. Afful E, Elliott B, Nayak MK, Phillips TW. Phosphine resistance in north American field populations of the lesser grain borer, Rhyzopertha dominica (Coleoptera: Bostrichidae). J Econ Entomol. 2018;111:463–9.

    Article  CAS  PubMed  Google Scholar 

  51. Cato A, Afful E, Nayak MK, Phillips TW. Evaluation of knockdown bioassay methods to assess phosphine resistance in the red flour beetle, Tribolium castaneum (Herbst) (Coleoptera: Tenebrionidae). Insects. 2019;10:140.

    Article  PubMed Central  Google Scholar 

  52. Nayak MK, Daglish GJ, Phillips TW. Managing resistance to chemical treatments in stored products pests. Stewart Postharvest Rev. 2015.

  53. Collins PJ. Resistance to grain protectants and fumigants in insect pests of stored products in Australia. In: Banks HJ, Wright EJ, Damcevski KA, editors. Stored Grain in Australia. Proceedings of the Australian Post-harvest Technical Conference, Canberra. Canberra: CSIRO, 26–29 May 1998; 1998. p. 55–7.

    Google Scholar 

  54. Vayias BJ, Athanassiou CG, Kavallieratos NG, Tsesmeli CD, Buchelos CT. Persistence and efficacy of two diatomaceous earth formulations and a mixture of diatomaceous earth with natural pyrethrum against Tribolium confusum Jacquelin du Val (Coleoptera: Tenebrionidae) on wheat and maize. Pest Manag Sci. 2006;62:456–64.

    Article  CAS  PubMed  Google Scholar 

  55. Yang Y, Wu Y, Chen S, Devine GJ, Denholm I, Jewess P, Moores GD. The involvement of microsomal oxidases in pyrethroid resistance in Helicoverpa armigera from Asia. Insect Biochem Mol Biol. 2004;34:763–73.

    Article  CAS  PubMed  Google Scholar 

  56. Nicolas J, Vincent P, Tralongo H, Lanteri P, Longeray R. Experimental design and modelization as tools for the development of insecticide combination of deltamethrin, organophosphorous compounds and piperonyl butoxide as stored grain protectant in temperate and tropical situations. In: FleuratLessard F, Ducom P, editors. Proceedings of the Fifth International Working Conference on Stored-product Protection. 9–14 September 1990, Bordeaux, France. Blanquefort Cedex: Imprimerie Medocaine; 1991. p. 589–98.

    Google Scholar 

  57. Yunta C, Grisales N, Nász S, Hemmings K, Pignatelli P, Voice M, Ranson H, Paine MJ. Pyriproxyfen is metabolized by P450s associated with pyrethroid resistance in an. Gambiae. Insect Biochem Mol Biol. 2016;78:50–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  58. Terriere LC, Yu SJ. Insect juvenile hormones: induction of detoxifying enzymes in the housefly and detoxification by housefly enzymes. Pestic Biochem Physiol. 1973;3:96–107.

    Article  CAS  Google Scholar 

  59. Giraudo M, Hilliou F, Fricaux T, Audatnt P, Feyereisen R, Le Goff G. Cytochrome P450s from the fall armyworm (Spodoptera frugiperda): responses to plant allelochemicals and pesticides. Insect Mol Biol. 2015;24:115–28.

    Article  CAS  PubMed  Google Scholar 

  60. Feyereisen R, Insect CYP. Genes and P450 enzymes. Insect Biochem Molec. 2012:236–316.

  61. Andersen JF, Walding JK, Evans PH, Bowers WS, Feyereisen R. Substrate specificity for the epoxidation of terpenoids and active site topology of house fly cytochrome P450 6A1. Chem Res Toxicol. 1997;10:156–64.

    Article  CAS  PubMed  Google Scholar 

  62. Karatolos N, Williamson MS, Denholm I, Gorman K, Ffrench-Constant RH, Bass C. Over-expression of a cytochrome P450 is associated with resistance to pyriproxyfen in the greenhouse whitefly Trialeurodes vaporariorum. PLoS One. 2012;7:e31077.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  63. Pfaffl MW, Horgan GW, Dempfle L. Relative expression software tool (REST©) for group-wise comparison and statistical analysis of relative expression results in real-time PCR. Nucleic Acids Res. 2002;30:e36.

    Article  PubMed  PubMed Central  Google Scholar 

  64. Grabherr MG, Haas BJ, Yassour M, Levin JZ, Thompson DA, Amit I, Adiconis X, Fan L, Raychowdhury R, Zeng Q, Chen Z, Mauceli E, Hacohen N, Gnirke A, Rhind N, di Palma F, Birren BW, Nusbaum C, Lindblad-Toh K, Friedman N, Regev A. Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nat Biotechnol. 2011;29:644–52.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  65. Jones P, Binns D, Chang HY, Fraser M, Li W, McAnulla C, McWilliam H, Maslen J, Mitchell A, Nuka G, Pesseat S, Quinn AF, Sangrador-Vegas A, Scheremetjew M, Yong SY, Lopez R, Hunter S. InterProScan 5: genome-scale protein function classification. Bioinformatics. 2014;30:1236–40.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  66. Camacho C, Coulouris G, Avagyan V, Ma N, Papadopoulos J, Bealer K, Madden TL. BLAST+: architecture and applications. BMC Bioinformatics. 2009;10:421.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  67. Suzek BE, Wang Y, Huang H, McGarvey PB, Wu CH, UniProt Consortium. UniRef clusters: a comprehensive and scalable alternative for improving sequence similarity searches. Bioinformatics. 2015;31:926–32.

    Article  CAS  PubMed  Google Scholar 

  68. Bray NL, Pimentel H, Melsted P, Pachter L. Near-optimal probabilistic RNA-seq quantification. Nat Biotechnol. 2016;34:525–7.

    Article  CAS  PubMed  Google Scholar 

  69. McCarthy DJ, Chen Y, Smyth GK. Differential expression analysis of multifactor RNA-Seq experiments with respect to biological variation. Nucleic Acids Res. 2012;40:4288–97.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  70. Reimand J, Arak T, Adler P, Kolberg L, Reisberg S, Peterson H, Vilo Z. G: profiler - a web server for functional interpretation of gene lists. Nucleic Acids Res. 2016;44:W83–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  71. Kim D, Langmead B, Salzberg SL. HISAT: a fast spliced aligner with low memory requirements. Nat Methods. 2015;12:357–60.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  72. Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, Marth G, Abecasis G. Durbin R; 1000 genome project data processing subgroup. The sequence alignment/map format and SAMtools. Bioinformatics. 2009;25:2078–9.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  73. Koboldt DC, Zhang Q, Larson DE, Shen D, MCLellan MD, Lin L, Miller CA, Mardis ER, Ding L, Wilson RK. VarScan 2: Somatic mutation and copy number alteration discovery in cancer by exome sequencing. Genome Res 2012; 22:568–576.

  74. Robinson JT, Thorvaldsdóttir H, Winckler W, Guttman M, Lander ES, Getz G, Mesirov JP. Integrative genomics viewer. Nat Biotechnol. 2011;29:24–6.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  75. Nelson DR. Cytochrome P450 and the individuality of species. Arch Biochem Biophys. 1999;369:1–10.

    Article  CAS  PubMed  Google Scholar 

  76. Katoh K, Standley DM. MAFFT multiple sequence alignment software version 7: improvements in performance and usability. Mol Biol Evol. 2013;30:772–80.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  77. Capella-Gutiérrez S, Silla-Martínez JM, Gabaldón T trimAl: A tool for automated alignment trimming in large-scale phylogenetic analyses Bioinformatics 2009; 25: 1972–1973.

  78. Stamatakis A. RAxML version 8: a tool for phylogenetic analysis and post-analysis of large phylogenies. Bioinformatics. 2014;30:1312–3.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  79. Stoer BC, Muller KF. TreeGraph 2: combining and visualizing evidence from different phylogenetic analyses. BMC Bioinformatics. 2010;11:7.

    Article  Google Scholar 

  80. He Z, Zhang H, Gao S, Lercher MJ, Chen WH, Hu S. Evolview v2: an online visualization and management tool for customized and annotated phylogenetic trees. Nucleic Acids Res. 2016;44:W236–41.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references


The authors would like to thank Dr. David Nelson for manually annotating the CYPs identified in this study.


This research has been co-financed by the European Union and Greek national funds through the Operational Program Competitiveness, Entrepreneurship and Innovation, under the call RESEARCH – CREATE – INNOVATE (project code: T1EDK-01491, “CropFeed”) to CGA and by the General Secretariat for Research and Technology of Greece (GSRT) grant “Plant-Up” (OPS5002803 - University of Crete) co-funded by GSRT and the EU Operational Programme “Competitiveness, Entrepreneurship and Innovation” (EPAnEK) within the National Strategic Reference Framework (NSRF 2014–2020), as well as a GSRT “Kripis” research program (to IMBB-FORTH) for the development of research institutes to JV.

Author information

Authors and Affiliations



CGA and JV designed the study. MR, PI contributed to the sequencing. GVB carried out the experiment procedure and MKS performed the analyses. MT performed the qPCR validation and MT and MR performed the analysis. CGA, JV, MR, PI, MKS, RJ and MKN wrote the manuscript and all authors read and made comments on the manuscript prior to submission.

Corresponding author

Correspondence to Maria K. Sakka.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Fig. S1.

Mean mortality (±SEM) of Rhyzopertha dominica after 7, 14 and 21 days of exposure for the susceptible (Lab-S) and resistant (Met-R) strain for all the combinations tested (control, 0.01 mg/kg, 0.03 mg/kg, 0.1 mg/kg, 0.3 mg/kg, PBO, 0.01 mg/kg + PBO, 0.03 mg/kg + PBO, 0.1 mg/kg + PBO, 0.3 mg/kg + PBO for susceptible and control, 1 mg/kg, 3 mg/kg, 10 mg/kg, 30 mg/kg, PBO, 1 mg/kg + PBO, 3 mg/kg + PBO, 10 mg/kg + PBO, 30 mg/kg + PBO for resistant).

Additional file 2: Fig. S2.

Principal components analysis of the transcript expression levels between the resistant (Met-R - blue) and susceptible (Lab-S - red) strains. The two strains are clearly different from each other, a prerequisite for downstream analyses.

Additional file 3: Table S1.

Corresponding information regarding sequence identity of CYPs.

Additional file 4: Table S2.

Details of all up- and down-regulated genes between the resistant Met-R and the susceptible Lab-S strains (padj < 0.001, log2|FC| > 2).

Additional file 5: Table S3.

Amino acid sequences used in the phylogenetic study of the R. dominica P450s, presented in Fig. 4.

Additional file 6: Table S4.

Primers used for the validation qPCRs.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Sakka, M.K., Riga, M., Ioannidis, P. et al. Transcriptomic analysis of s-methoprene resistance in the lesser grain borer, Rhyzopertha dominica, and evaluation of piperonyl butoxide as a resistance breaker. BMC Genomics 22, 65 (2021).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: