==== Front Sci Rep Sci Rep Scientific Reports 2045-2322 Nature Publishing Group UK London 33319828 78974 10.1038/s41598-020-78974-z Article Selection and validation of experimental condition-specific reference genes for qRT-PCR in Metopolophium dirhodum (Walker) (Hemiptera: Aphididae) Li Xinan 12 Gong Peipan 1 Wang Bingting 3 Wang Chao 1 Li Mengyi 1 Zhang Yunhui 1 Li Xiangrui 1 Gao Haifeng 4 Ju Jiansong jujiansong@126.com 3 Zhu Xun zhuxun@caas.cn 1 1 grid.410727.70000 0001 0526 1937State Key Laboratory for Biology of Plant Diseases and Insect Pests, Institute of Plant Protection, Chinese Academy of Agricultural Sciences, Yuanmingyuan West Road No.2, Haidian District, Beijing, 100193 China 2 grid.503006.00000 0004 1761 7808School of Resource and Environmental Sciences, Henan Institute of Science and Technology, Eastern HuaLan Avenue, Xinxiang, 453003 Henan China 3 grid.256884.50000 0004 0605 1239College of Life Science, Hebei Normal University, Road Nan er huan dong No.20, Shijiazhuang, 050024 Hebei China 4 grid.418524.e0000 0004 0369 6250Institute of Plant Protection, Xinjiang Academy of Agricultural Sciences/Key Laboratory of Integrated Pest Management on Crop in Northwestern Oasis, Ministry of Agriculture and Rural Affairs, Ürümqi, 830091 China 15 12 2020 15 12 2020 2020 10 219519 2 2020 5 11 2020 © The Author(s) 2020Open 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 http://creativecommons.org/licenses/by/4.0/.Metopolophium dirhodum (Walker) (Hemiptera: Aphididae) is one of the most common aphid pests of winter cereals. To facilitate accurate gene expression analyses with qRT-PCR assays, the expression stability of candidate reference genes under specific experimental conditions must be verified before they can be used to normalize target gene expression levels. In this study, 10 candidate reference genes in M. dirhodum were analyzed by qRT-PCR under various experimental conditions. Their expression stability was evaluated with delta Ct, BestKeeper, geNorm, and NormFinder methods, and the final stability ranking was determined with RefFinder. The results indicate that the most appropriate sets of internal controls were SDHB and RPL8 across geographic population; RPL8, Actin, and GAPDH across developmental stage; SDHB and NADH across body part; RPL8 and Actin across wing dimorphism and temperature; RPL4 and EF1A across starvation stress; AK and RPL4 across insecticide treatments; RPL8 and NADH across antibiotic treatments; RPL8, RPL4, Actin, and NADH across all samples. The results of this study provide useful insights for establishing a standardized qRT-PCR procedure for M. dirhodum and may be relevant for identifying appropriate reference genes for molecular analyses of related insects. Subject terms Reverse transcription polymerase chain reactionEntomologyChina Agriculture Research SystemCARS-3National Key Research and Development Program of China2017YFD02017002018YFD0200500Zhang Yunhui Li Xiangrui issue-copyright-statement© The Author(s) 2020 ==== Body Introduction The quantitative analysis of target gene expression is an essential part of most molecular studies. Quantitative real-time PCR (qRT-PCR) is a powerful tool for quantifying gene expression, combining improvements in both sensitivity and specificity with efficient techniques for signal detection. It is useful for the quantitative data analysis required for research related to molecular medicine, biotechnology, microbiology, and diagnostics and has become the preferred method for quantifying mRNA1. Nevertheless, gene expression analyses are affected by many factors such as the quality of RNA samples, the efficiency of reverse transcription, and PCR efficiency2,3. For accurate comparisons of expression levels, the expression data of the genes of interest are normalized against the expression data for a reference gene4. Moreover, the reference gene compensates for the above-mentioned limitations5. Because housekeeping genes are related to ubiquitous and basic cellular functions, they are considered to be constitutively expressed under diverse conditions6. Housekeeping genes, including those encoding actin, glyceraldehyde-3-phosphate dehydrogenase, ribosomal protein, 18S ribosomal RNA, elongation factor 1α and heat shock proteins, have been extensively used as endogenous controls for normalizing real-time PCR data7–11. However, several studies have indicated that the expression levels of the reference genes vary under diverse conditions12–14. In fact, no single reference gene is appropriate for all experimental conditions. Therefore, evaluating and validating the stability of reference genes under different experimental conditions is critical. There have recently been several reports regarding reference genes for molecular research on insects, including bumblebee, Harmonia axyridis, Propylea japonica, Aphis craccivora Koch, Henosepilachna vigintioctomaculata, Chilo suppressalis, Galeruca daurica, Liriomyza trifolii, Coccinella septempunctata, Phenacoccus solenopsis, Lipaphis erysimi, Myzus persicae, Acyrthosiphon pisum, and Megoura viciae11,14–27. Metopolophium dirhodum (Walker) (Hemiptera: Aphididae) is one of the most major aphid pests affecting winter wheat and other cereals worldwide28–31. Additionally, M. dirhodum, which was first detected in the 1970s, originated in the Holarctic and was subsequently introduced to South America and other regions32,33. The M. dirhodum nymphs and adults damage cereals by directly feeding on plants, which may result in grain yield losses of 27–30%34. Moreover, they damage crops by transmitting several viruses, especially the barley yellow dwarf virus35. This aphid has most often been detected in semi-arid regions in South America, South Africa, Australia, and New Zealand, where it damages cereals, including wheat, barley, rye, and oat. A previous study revealed that M. dirhodum is the most abundant aphid species on cereals in the continental climate of central Europe33. With the technical advances occurring in the post-genomic era, researchers may soon have additional options for studying M. dirhodum at the molecular level, which may contribute to the development of improved control measures. Thus, identifying suitable reference genes is important for analyzing the expression of functional genes and for evaluating the efficiency of target gene silencing via RNA interference. The objective of this study was to identify and evaluate a suite of experimental condition-specific reference genes to normalize target gene expression in M. dirhodum. Specifically, we analyzed the following 10 candidate genes: Actin, glyceraldehyde-3-phosphate dehydrogenase (GAPDH), NADH dehydrogenase (NADH), arginine kinase (AK), succinate dehydrogenase B (SDHB), ribosomal protein L8 (RPL8), 18S ribosomal RNA (18S), elongation factor 1α (EF1A), ribosomal protein L4 (RPL4), and heat shock protein 68 (HSP68). The effects of the following factors on reference gene expression were evaluated: geographic population, developmental stage, body part, wing dimorphism, temperature, starvation stress, and exposure to an insecticide or antibiotic. The results indicate that the best reference genes for analyzing M. dirhodum gene expression vary among conditions. Results Expression levels of candidate reference genes To evaluate the expression profiles of the selected candidate genes in all M. dirhodum sample sets, mRNA levels were measured for all genes. The gene expression levels varied considerably between Ct values of 12.70 (18S) and 30.88 (GAPDH) (Fig. 1). Of the 10 analyzed genes, the highest and lowest expression levels were detected for 18S (mean Ct value of 14.27) and GAPDH (mean Ct value of 28.90), respectively. The least variable expression among all samples was observed for Actin (mean Ct value ± SD of 26.79 ± 0.42) and RPL8 (21.10 ± 0.35). In contrast, HSP68 (24.82 ± 1.86) exhibited the most variable expression in all the tested samples.Figure 1 Candidate reference gene expression levels. Candidate reference gene expression levels in the whole M. dirhodum sample set are expressed in terms of the threshold cycle number (Ct value). Data are presented as whisker box plots. The box represents the 25th–75th percentiles, the median is indicated by a bar across the box, and the whiskers on each box represent the minimum and maximum values. Analysis of gene expression stability Geographic populations The delta Ct method and the BestKeeper, NormFinder, and geNorm algorithms were used to assess the stability of the candidate reference gene expression levels. The rank order (most to least stable expression) was highly consistent among the four methods. Specifically, SDHB, RPL4, and RPL8 were identified as the most stable genes, whereas HSP68 and GAPDH were the least stable genes (Table 1). The RefFinder results for the geographic populations revealed a rank order (most to least stable expression) of SDHB, RPL8, RPL4, NADH, AK, 18S, Actin, EF1A, GAPDH, and HSP68 (Fig. 2). On the basis of the GeNorm analysis, all pairwise variation values were below the 0.15 cut-off value, except for V5/6 (Fig. 3). Moreover, the RefFinder analysis indicated SDHB and RPL8 are required for the normalization of target gene expression levels in different geographic populations.Table 1 Rank order of the M. dirhodum candidate reference genes under various experimental conditions. Experimental conditions Rank Delta CT BestKeeper NormFinder GeNorm Gene name Standard deviation Gene name Standard deviation Gene name Stability value Gene name Stability value Geographic populations 1 SDHB 0.73 SDHB 0.04 RPL4 0.292 SDHB/RPL8 0.123 2 RPL4 0.74 RPL8 0.12 SDHB 0.311 3 NADH 0.78 NADH 0.13 RPL8 0.439 NADH 0.129 4 RPL8 0.78 RPL4 0.20 NADH 0.474 RPL4 0.225 5 AK 0.85 AK 0.29 18S 0.537 AK 0.274 6 18S 0.90 18S 0.62 AK 0.593 18S 0.507 7 Actin 0.96 Actin 0.68 Actin 0.707 Actin 0.640 8 EF1A 1.06 EF1A 0.75 EF1A 0.862 EF1A 0.727 9 GAPDH 1.38 HSP68 0.84 GAPDH 1.310 GAPDH 0.843 10 HSP68 1.44 GAPDH 1.08 HSP68 1.366 HSP68 0.962 Development-al stages 1 GAPDH 1.03 RPL8 0.61 GAPDH 0.149 Actin/RPL8 0.461 2 Actin 1.04 RPL4 0.64 Actin 0.231 3 RPL8 1.08 Actin 0.85 RPL8 0.410 GAPDH 0.504 4 NADH 1.10 SDHB 0.89 NADH 0.458 NADH 0.564 5 RPL4 1.31 NADH 0.94 RPL4 0.934 RPL4 0.653 6 18S 1.55 GAPDH 0.96 AK 1.199 EF1A 0.840 7 AK 1.55 AK 1.40 18S 1.266 18S 0.949 8 EF1A 1.56 18S 1.47 EF1A 1.301 AK 1.123 9 SDHB 1.75 EF1A 1.50 SDHB 1.503 SDHB 1.251 10 HSP68 1.97 HSP68 1.57 HSP68 1.763 HSP68 1.395 Body parts 1 SDHB 0.72 GAPDH 0.24 NADH 0.043 NADH/SDHB 0.085 2 NADH 0.74 18S 0.28 SDHB 0.043 3 18S 0.77 EF1A 0.32 18S 0.168 Actin 0.270 4 Actin 0.81 SDHB 0.36 Actin 0.447 18S 0.314 5 EF1A 0.88 NADH 0.42 EF1A 0.568 AK 0.377 6 AK 0.95 Actin 0.50 GAPDH 0.687 EF1A 0.554 7 GAPDH 0.98 RPL8 0.61 AK 0.711 GAPDH 0.645 8 RPL8 1.15 AK 0.66 RPL8 1.056 RPL8 0.757 9 RPL4 1.25 RPL4 0.71 RPL4 1.162 RPL4 0.832 10 HSP68 1.59 HSP68 1.23 HSP68 1.557 HSP68 0.983 Wing dimorphism 1 Actin 0.60 Actin 0.06 RPL8 0.027 RPL8/EF1A 0.053 2 RPL4 0.60 RPL4 0.14 EF1A 0.027 3 RPL8 0.62 HSP68 0.19 RPL4 0.039 RPL4 0.087 4 EF1A 0.64 RPL8 0.19 Actin 0.094 Actin 0.217 5 HSP68 0.64 NADH 0.23 HSP68 0.332 HSP68 0.309 6 NADH 0.67 EF1A 0.23 NADH 0.406 NADH 0.344 7 SDHB 0.94 SDHB 0.51 SDHB 0.889 SDHB 0.458 8 AK 1.01 AK 0.57 AK 0.982 AK 0.523 9 GAPDH 1.13 GAPDH 0.74 GAPDH 1.035 GAPDH 0.680 10 18S 1.41 18S 0.97 18S 1.401 18S 0.827 Temperatures 1 Actin 0.72 RPL8 0.10 RPL4 0.032 Actin/NADH 0.206 2 RPL8 0.74 RPL4 0.16 RPL8 0.064 3 NADH 0.75 SDHB 0.29 Actin 0.141 RPL8 0.280 4 RPL4 0.78 Actin 0.30 EF1A 0.266 RPL4 0.310 5 EF1A 0.81 EF1A 0.30 NADH 0.347 EF1A 0.341 6 SDHB 0.83 NADH 0.32 SDHB 0.397 SDHB 0.386 7 AK 0.94 AK 0.46 AK 0.502 AK 0.448 8 GAPDH 1.00 GAPDH 0.60 GAPDH 0.626 GAPDH 0.526 9 18S 1.11 18S 0.69 18S 0.915 18S 0.595 10 HSP68 2.92 HSP68 2.18 HSP68 2.885 HSP68 1.059 Starvation-stress 1 RPL4 1.03 18S 0.06 EF1A 0.026 NADH/AK 0.050 2 EF1A 1.03 Actin 0.33 RPL4 0.026 3 RPL8 1.10 GAPDH 0.46 RPL8 0.484 SDHB 0.175 4 AK 1.23 RPL8 0.78 AK 0.706 RPL4 0.599 5 NADH 1.26 EF1A 1.02 NADH 0.771 EF1A 0.687 6 GAPDH 1.31 RPL4 1.05 SDHB 1.034 RPL8 0.790 7 SDHB 1.40 AK 1.70 GAPDH 1.066 GAPDH 0.930 8 Actin 1.43 NADH 1.74 Actin 1.280 Actin 1.017 9 18S 1.77 SDHB 1.89 18S 1.727 18S 1.125 10 HSP68 2.55 HSP68 2.80 HSP68 2.527 HSP68 1.409 Insecticide-stress 1 RPL4 0.32 HSP68 0.12 AK 0.129 Actin/AK 0.028 2 AK 0.32 SDHB 0.22 RPL4 0.135 3 Actin 0.33 RPL8 0.22 NADH 0.154 RPL8 0.080 4 RPL8 0.33 RPL4 0.23 Actin 0.167 RPL4 0.102 5 NADH 0.37 Actin 0.29 GAPDH 0.192 HSP68 0.151 6 GAPDH 0.39 NADH 0.29 RPL8 0.208 NADH 0.205 7 HSP68 0.44 AK 0.31 SDHB 0.384 SDHB 0.245 8 SDHB 0.47 GAPDH 0.50 HSP68 0.388 GAPDH 0.281 9 18S 0.56 18S 0.69 18S 0.478 18S 0.347 10 EF1A 0.76 EF1A 0.77 EF1A 0.731 EF1A 0.431 Antibiotic-stress 1 RPL8 0.54 SDHB 0.03 NADH 0.024 GAPDH/18S 0.013 2 RPL4 0.54 Actin 0.15 RPL8 0.086 3 AK 0.54 NADH 0.18 Actin 0.087 AK 0.060 4 18S 0.58 RPL8 0.47 RPL4 0.350 RPL4 0.071 5 GAPDH 0.59 RPL4 0.59 SDHB 0.371 EF1A 0.112 6 NADH 0.63 AK 0.61 AK 0.383 RPL8 0.163 7 Actin 0.65 18S 0.67 18S 0.484 NADH 0.297 8 EF1A 0.69 GAPDH 0.68 GAPDH 0.501 Actin 0.370 9 SDHB 0.76 EF1A 0.76 EF1A 0.646 SDHB 0.439 10 HSP68 1.99 HSP68 0.95 HSP68 1.987 HSP68 0.749 All above conditions 1 RPL8 1.01 Actin 0.54 RPL8 0.401 RPL8/RPL4 0.421 2 RPL4 1.03 RPL8 0.54 RPL4 0.497 3 NADH 1.09 RPL4 0.82 Actin 0.543 EF1A 0.674 4 Actin 1.10 18S 0.96 NADH 0.624 NADH 0.747 5 EF1A 1.15 SDHB 1.00 SDHB 0.723 GAPDH 0.786 6 GAPDH 1.16 EF1A 1.01 EF1A 0.724 Actin 0.827 7 SDHB 1.17 GAPDH 1.15 GAPDH 0.752 SDHB 0.868 8 AK 1.44 NADH 1.16 AK 1.159 AK 0.955 9 18S 1.51 HSP68 1.46 18S 1.230 18S 1.061 10 HSP68 2.16 AK 1.56 HSP68 2.019 HSP68 1.281 Figure 2 Stability of candidate reference gene expression levels in response to various treatments and conditions. In a RefFinder analysis, decreasing Geomean values correspond to increasing gene expression stability. The Geomean values for the following M. dirhodum samples are presented: adult samples from different geographic populations (Geographic population), samples for all developmental stages (Developmental stages), samples for different body parts of wingless adults (Body part), samples for winged and wingless adults (Wing dimorphism), adult samples exposed to different temperatures (Temperature-stress), fed and unfed adult samples (Starvation-stress), adult samples treated with different insecticides (Insecticide-stress), adult samples treated with antibiotic (Antibiotic-stress), and all samples for all treatments (All conditions). The candidate reference genes are as follows: Actin, Actin; GAPDH, glyceraldehyde-3-phosphate dehydrogenase; NADH, NADH dehydrogenase; AK, arginine kinase; SDHB, succinate dehydrogenase B; RPL8, ribosomal protein L18; RPL4, ribosomal protein L4; HSP68, heat shock protein 68; 18S, 18S ribosomal RNA; and EF1A, elongation factor 1α. Figure 3 Determination of the optimal number of reference genes for accurate normalization calculated by geNorm. The Vn/n+1 value indicates the pairwise variation (Y axis) between two sequential normalization factors and determines the optimal number of reference genes required for an accurate data normalization. A value below 0.15 indicates that an additional reference gene will not significantly improve the normalization. Developmental stage The delta Ct and NormFinder analyses identified GAPDH and Actin as the most stable genes. In contrast, the most stable genes were RPL8 and RPL4 according to BestKeeper and Actin and RPL8 according to GeNorm. Regardless of the method, HSP68 was identified as the least stable gene (Table 1). According to the RefFinder results, the rank order (most to least stable expression) for the developmental stages was RPL8, Actin, GAPDH, RPL4, NADH, 18S, AK, SDHB, EF1A, and HSP68 (Fig. 2). The GeNorm analysis revealed that the values for V3/4 were less than the proposed 0.15 cut-off (Fig. 3). The RefFinder analysis indicated RPL8, Actin, and GAPDH are required for normalizing target gene expression levels for the different M. dirhodum developmental stages. Body part The gene expression stability rank order determined with BestKeeper differed from that obtained with the other three methods (Table 1). The BestKeeper algorithm identified GAPDH and 18S as the most stable genes. In contrast, the delta Ct method, NormFinder, and GeNorm identified NADH and SDHB as the most stable genes. All four analyses indicated RPL4 and HSP68 were the least stable genes. The RefFinder results for the different body parts revealed a rank order (most to least stable expression) of SDHB, NADH, 18S, Actin, GAPDH, EF1A, AK, RPL8, RPL4, and HSP68 (Fig. 2). On the basis of the GeNorm analysis, all pairwise variation values were below the 0.15 cut-off value, except for V9/10 (Fig. 3). The RefFinder analysis indicated SDHB and NADH are required for normalizing target gene expression levels in various M. dirhodum body parts. Wing dimorphism The delta Ct and BestKeeper analyses identified Actin and RPL4 as the most stable genes, whereas both NormFinder and GeNorm identified RPL8 and EF1A as the most stable genes. All four analyses indicated that 18S, GAPDH, AK, and SDHB were the least stable genes (Table 1). The RefFinder data for the wing dimorphism revealed a rank order (most to least stable expression) of RPL8, Actin, RPL4, EF1A, HSP68, NADH, SDHB, AK, GAPDH, and 18S (Fig. 2). On the basis of the GeNorm analysis, all pairwise variation values were below the 0.15 cut-off value (Fig. 3). According to RefFinder, RPL8 and Actin are required for normalizing target gene expression levels in wing-dimorphic insects. Temperature-induced stress The delta Ct method identified Actin and RPL8 as the most stable genes. Both BestKeeper and NormFinder identified RPL8 and RPL4 as the most stable genes, whereas GeNorm identified Actin and NADH as the most stable genes. All four analyses indicated HSP68, 18S, GAPDH, and AK were the least stable genes (Table 1). The RefFinder data for the different temperatures revealed a rank order (most to least stable expression) of RPL8, Actin, RPL4, NADH, EF1A, SDHB, AK, GAPDH, 18S, and HSP68 (Fig. 2). On the basis of the GeNorm analysis, all pairwise variation values were below the 0.15 cut-off value, except for V9/10 (Fig. 3). The RefFinder analysis indicated RPL8 and Actin are required for normalizing target gene expression levels in M. dirhodum exposed to different temperatures. Starvation-induced stress The delta Ct method and the NormFinder algorithm identified EF1A and RPL4 as the most stable genes and Actin, 18S, and HSP68 as the least stable genes (Table 1). However, BestKeeper identified 18S and Actin as the most stable genes and SDHB and HSP68 as the least stable genes (Table 1). The GeNorm algorithm identified NADH and AK as the most stable genes and Actin, 18S, and HSP68 as the least stable genes (Table 1). The RefFinder results for the starvation treatment revealed a rank order (most to least stable expression) of RPL4, EF1A, AK, NADH, RPL8, 18S, GAPDH, Actin, SDHB, and HSP68 (Fig. 2). The GeNorm analysis indicated that the pairwise variation values for V2/3 were less than the proposed 0.15 cut-off (Fig. 3). The RefFinder analysis indicated RPL4 and EF1A are required for normalizing target gene expression levels in starvation-stressed M. dirhodum. Insecticide-induced stress The delta Ct and NormFinder data revealed AK and RPL4 as the most stable genes, whereas the BestKeeper results identified HSP68 and SDHB as the most stable genes. In contrast, Actin and AK were the most stable genes according to GeNorm. All four analyses identified 18S and EF1A as the least stable genes (Table 1). The RefFinder data for the insecticide treatment revealed a rank order (most to least stable expression) of AK, RPL4, Actin, RPL8, HSP68, NADH, SDHB, GAPDH, 18S, and EF1A (Fig. 2). Based on the GeNorm analysis, all the pairwise variation values were below 0.15 cut-off value (Fig. 3). Thus, AK and RPL4 are required for normalizing target gene expression levels in insecticide-treated M. dirhodum. Antibiotic-induced stress The delta Ct method identified RPL8 and RPL4 as the most stable genes. The BestKeeper algorithm identified SDHB and Actin as the most stable genes, whereas NormFinder indicated NADH and RPL8 were the most stable genes. The GeNorm algorithm identified GAPDH and 18S as the most stable genes. All four analyses identified EF1A, SDHB, and HSP68 as the least stable genes (Table 1). The RefFinder data for the antibiotic treatment revealed a rank order (most to least stable expression) of RPL8, NADH, RPL4, 18S, GAPDH, AK, Actin, SDHB, EF1A, and HSP68 (Fig. 2). According to the GeNorm analysis, all pairwise variation values were less than the proposed 0.15 cut-off, except for V9/10 (Fig. 3). The RefFinder analysis suggested RPL8 and NADH are required for normalizing the target gene expression levels in antibiotic-treated M. dirhodum. Overall ranking of M. dirhodum candidate reference genes An examination of the candidate reference gene expression stability for all treatments and conditions with the four methods used in this study produced similar rank orders, with RPL4 and RPL8 identified as the most stable genes and AK, 18S, and HSP68 revealed as the least stable genes (Table 1). The RefFinder results for all treatments and conditions revealed a rank order (most to least stable expression) of RPL8, RPL4, Actin, NADH, EF1A, SDHB, GAPDH, 18S, AK, and HSP68 (Fig. 2). The GeNorm analysis indicated that the pairwise variation values for V4/5 were less than the proposed 0.15 cut-off (Fig. 3). Thus, an analysis of all treatments and conditions suggested that RPL8, RPL4, Actin, and NADH are suitable internal reference genes for normalizing target gene expression levels in M. dirhodum. Discussion There are several reports describing the application of qRT-PCR assays to clarify the gene expression levels associated with diverse biological processes36–39. Reference genes used for molecular investigations can influence the accuracy of target gene expression levels6,40–42. Therefore, a stable reference gene is an important prerequisite for gene expression investigations. Housekeeping genes, which are constitutively expressed to maintain basic cellular functions, have traditionally been used as internal reference controls6,10,11. However, there is no universal reference gene that is stably expressed in all cell and tissue types under different experimental conditions10,11,43–47. Therefore, every stable reference gene used to normalize gene expression data should be evaluated under each experimental condition43,48. In this study, qRT-PCR was used to evaluate the expression-level stability of 10 candidate reference genes in M. dirhodum across specific conditions. The best reference genes varied among conditions. Specifically, RPL8 (mean Ct value ± SD, 21.10 ± 0.35) and Actin (26.79 ± 0.42) had the least variable expression levels, whereas HSP68 (24.82 ± 1.86) produced the most variable expression levels among the examined candidate reference genes (Fig. 1). Similarly, RPL8, RPL4, and Actin were the most stable reference genes, whereas HSP68 and 18S were the least stable reference genes under most conditions (Fig. 2). Ribosomal proteins (RPs), which are the principal components of ribosomes, are one of the most highly conserved proteins in all life forms. Earlier research proved that RP-encoding genes are among the most stably expressed reference genes, and have been widely used to normalize gene expression levels in insect molecular investigations during the past 10 years49. For example, in Bradysia odoriphaga50, RPS15 was the most stably expressed gene in response to various temperature treatments. However, another study indicated that the expression levels of RP-encoding genes may vary under some conditions49. Moreover, RPS20 was detected as the least stably expressed gene for analyzing Plutella xylostella geographic populations as well as the effects of the temperature, photoperiod, and insecticides10. Consistent with these earlier findings, we identified RPL8 as the most stable gene in M. dirhodum across various conditions (except for analyses of different body parts, starvation stress, and insecticide treatments) (Fig. 2). Additionally, RPL4 was detected as the most stable gene in response to starvation and insecticide treatments, but was also almost the least stable gene during analyses of various M. dirhodum body parts (Fig. 2). Actin, which encodes a major structural protein, is important for cell secretion, motility, cytoplasm flow, and cytoskeleton maintenance. Moreover, Actin is expressed at various levels in many cell types, and is considered the ideal reference gene for qRT-PCR, which may explain its frequent use15,26. For example, it has been used to study the effects of diet on B. odoriphaga gene expression50 and for investigating M. persicae gene expression in different tissues and in response to the temperature, photoperiod, and wing dimorphism26. However, in Helicoverpa armigera, Actin was revealed to be the least stable reference gene following temperature and photoperiod treatments51. In our study, Actin was identified as one of the most stable reference genes for analyzing developmental stages, temperature effects, and wing dimorphism (Fig. 2). The GAPDH gene has been commonly used as a reference gene in the studies of gene expression7,52,53. However, unstable GAPDH expression has been detected in Tetranychus cinnabarinus developmental stages54, in the labial glands and fat bodies of Bombus terrestris and Bombus lucorum55, and in various Sogatella furcifera body parts56. In the current study, GAPDH was revealed as a stably expressed candidate reference gene for analyses of developmental stages (Fig. 2). These results imply that the mechanism underlying the expression stability of endogenous reference genes is complex. Furthermore, the stability of potential reference genes in different biological samples should be tested prior to their use. The protein encoded by EF1A affects translation by catalyzing the GTP-dependent binding of aminoacyl-tRNA to the acceptor site of the ribosome. The EF1A gene was recently used as a reference gene in multiple insect gene expression studies55,57,58. Our results suggest that EF1A is an appropriate reference gene only for analyzing the effects of starvation stress on M. dirhodum gene expression (Table 2).Table 2 Recommended reference genes for M. dirhodum under various experimental conditions. Conditions Reference gene Conditions Reference gene Population SDHB, RPL8 Temperature RPL8, Actin Development stage RPL8, Actin, GAPDH Starvation RPL4, EF1A Body part SDHB, NADH Insecticide AK, RPL4 Wing dimorphism RPL8, Actin Antibiotic RPL8, NADH All conditions RPL8, RPL4, Actin, NADH The AK gene encodes the phosphagen kinase in invertebrates, and it has rarely been used as a reference gene59. An earlier study revealed that AK is the most stably expressed gene in the B. terrestris labial gland and fat body60. In this study, AK was identified as the most stable gene following insecticide treatments (Fig. 2). In A. pisum, SDHB and NADH are reportedly the most stable housekeeping genes in developmental stages and in response to various temperatures11. However, we determined that SDHB and NADH are the most stable housekeeping genes only during examinations of different M. dirhodum body parts (Fig. 2). These results further suggest that reference gene expression stability is influenced by the experimental conditions. The 18S rRNA gene is considered to be an ideal reference control because of its relatively stable expression levels61. Accordingly, it has been applied in previous studies involving Lucilia cuprina62, Rhodnius prolixus63,64, and Delphacodes kuscheli65. However, in this study, 18S was revealed as one of the least stable genes in almost all sample sets, implying it is an inappropriate reference gene for M. dirhodum (Fig. 2). This observation is consistent with the results of previous studies that indicated 18S rRNA is not a stable reference gene in Bactrocera dorsalis and Nilaparvata lugens under specific experimental conditions66. It is transcribed by a separate RNA polymerase, which may explain why rRNA is not a suitable reference control67. Moreover, the utility of 18S for normalizing target gene expression levels in a qRT-PCR assay is limited by the potential imbalance between rRNA and mRNA fractions among samples61. The HSP68 gene, which belongs to the HSP70 family, encodes a highly conserved chaperone involved in protein assembly, folding, and transport as well as in antigen processing and presentation. The expression of genes encoding HSPs can be affected by high temperatures or other stresses (e.g., due to chemicals)68. In the current study, HSP68 was the least stable gene for all conditions (Fig. 2). In a previous study on Coleomegilla maculata, HSP70 was identified as the most stably expressed gene for sexes, but was the least stably expressed gene for analyses of different tissues, and dsRNA exposure44. It is becoming common for researchers to use multiple reference genes to normalize target gene expression levels in diverse studies because a single gene is usually insufficient for analyzing gene expression69. An earlier investigation indicated that too many or too few reference genes may adversely affect the robustness of data normalizations70. However, the simultaneous application of multiple reference genes in a given experiment may decrease the probability of biased normalizations. The optimal number of reference genes under specific experimental conditions can be determined with the geNorm algorithm, which calculates the pairwise variation Vn/n+1 based on the normalization factors NFn and NFn+1, with n ≥ 2. If Vn/n+1 is below 0.15, n is the optimal number of reference genes. The results of this study indicate that the most appropriate number of reference genes varies under diverse experimental conditions (Fig. 3). This implies that the stability of reference genes must be evaluated before every qRT-PCR experiment. Conclusions To the best of our knowledge, this study is the first to evaluate and validate experimental condition-specific candidate reference genes for M. dirhodum gene expression analyses. We identified reference genes applicable for elucidating functional gene expression profiles. In this study, we examined 10 candidate reference genes under diverse conditions. Notably, the stability of candidate gene expression levels in M. dirhodum varies depending on the experimental conditions. Moreover, we identified internal reference genes suitable for normalizing and quantifying gene expression in M. dirhodum (Table 2). Our findings may be useful for establishing a more accurate and reliable method for normalizing M. dirhodum qRT-PCR data. They may also provide the basis for future investigations on RNA interference and gene transcription in M. dirhodum and other insect pests. Materials and methods Insects Our original M. dirhodum colony was collected in Yinchuan (Ningxia), China (38° 48′ 54.78″ N, 106° 30′ 27.93″ E) in 2018. Other colonies were collected in Langfang (Hebei), China (39° 8′ 9.8″ N, 116° 10′ 4.05″ E) and Guiyang (Guizhou), China (26° 0′ 34.08″ N, 106° 35′ 4.35″ E) in 2018. The alive adults were collected in wheat leaves of different plants of these geographic locations and were taken back to the lab to establish population. All the wheat aphid populations were reared on Lunxuan 987 wheat seedlings in a thermostatic chamber maintained at 20 ± 2 °C and 60% relative humidity, with a 16-h light:8-h dark cycle. Treatments Geographic population Insects collected in Yinchuan (Ningxia), Langfang (Hebei), and Guiyang (Guizhou) in 2018 were examined to assess the effects of geography on gene expression. These three locations are separated by more than 1000 km. For each geographic population, three samples of 20 adults were selected, flash frozen in liquid nitrogen, and stored at − 80 °C until total RNA extraction. Developmental stage Three M. dirhodum samples of about 30 first-instar nymphs, 30 second-instar nymphs, 20 third-instar nymphs, 20 fourth-instar nymphs, and 20 adults were collected, flash frozen in liquid nitrogen, and stored at − 80 °C until total RNA extraction. Body part We used a dissection needle and a tweezer to separate the head, thorax, and abdomen from wingless M. dirhodum adults. These body parts as well as whole adult bodies were stored as described earlier. Wing dimorphism Three samples of 20 winged and wingless M. dirhodum adults were collected, flash frozen in liquid nitrogen, and stored at − 80 °C until total RNA extraction. Temperature-induced stress Potted wheat seedlings infested with M. dirhodum were divided into five groups for a 24-h exposure to one of the following five temperatures: 4, 10, 15, 20, and 25 °C. For each temperature, three samples of 20 adults were collected, flash frozen in liquid nitrogen, and stored at − 80 °C until total RNA extraction. None of the temperature treatments were lethal to the aphids. Starvation-induced stress Adult aphids were placed on moistened filter paper in a Petri dish (9 cm diameter) with no food for a 32-h incubation in a thermostatic chamber (20 ± 2 °C and 60% relative humidity, with a 16-h light:8-h dark cycle). The control (satiated) group comprised aphids able to feed on wheat seedlings in the same conditions. For the control and treatment groups, three samples of 20 adults were collected, flash frozen in liquid nitrogen, and stored at − 80 °C until total RNA extraction. The mortality rate among the starved aphids was approximately 10%. Insecticide-induced stress The effects of insecticides on the stability of candidate reference genes were assessed in M. dirhodum subjected to one of the following three insecticide treatments: imidacloprid (9.87 mg/L), thiamethoxam (122.00 mg/L), and beta-cypermethrin (17.28 mg/L). These concentrations were selected because a bioassay indicated they are 30% to the mortality of the population (LC30) (Table S1). Aphids were treated with the insecticides via the leaf dip method71. The 1% insecticide stock solutions prepared in acetone were serially diluted with water (containing 0.1% Tween-80) to produce five concentrations. Water (containing 0.1% Tween-80) was used as a control solution. Wheat leaves with M. dirhodum were immersed in the prepared solutions for 3–5 s and then placed on moistened filter paper in a Petri dish (9 cm diameter). The samples were incubated for 24 h at 20 ± 2 °C and 60% relative humidity, with a 16-h light:8-h dark cycle. For each concentration, the mortality rate based on three replicates of 30 aphids was calculated. Additionally, for the control and treatment groups, three samples of 20 adults were collected, flash frozen in liquid nitrogen, and stored at − 80 °C until total RNA extraction. Antibiotic-induced stress The M. dirhodum adults were fed a 30% sucrose solution containing 50 µg/mL rifampicin or an antibiotic-free sucrose solution (control) (25 aphids per feeder) for 48 h72. For the control and treatment groups, three samples of 20 adults were collected, flash frozen in liquid nitrogen, and stored at − 80 °C until total RNA extraction. Total RNA extraction and cDNA synthesis Total RNA was extracted with Trizol according to the protocol for the TRNzol Universal Reagent (Tiangen, Beijing, China). The ratio of the absorbance at 260 and 280 nm was 1.981–2.121, indicating the extracted RNA was pure. Next, 1 μg RNA was used as the template to synthesize first-strand cDNA with Oligo dT primers using the FastKing gDNA Dispelling RT SuperMix (Tiangen) following the manufacturer-recommended protocol. The synthesized cDNA was stored at − 20 °C. Primer design and quantitative real-time PCR A qRT-PCR assay was completed with the Talent qPCR PreMix (SYBR Green; Tiangen) and the CFX Connect Real-Time system (Bio-Rad, Hercules, CA, USA). Details regarding the primers for EF1A and 18S (Table 3) have been published by NCBI. Primers for the other target genes were designed based on our unpublished RNA sequencing data for M. dirhodum. The cDNA of each sample was prepared as a 50 ng/μL working solution. The qRT-PCR was completed in a 25-μL reaction volume comprising 12.5 μL 2 × Talent qPCR PreMix, 1 μL forward primer (100 μM), 1 μL reverse primer (100 μM), 1 μL cDNA working solution, and 9.5 μL RNase-Free ddH20. The PCR program was as follows: 95 °C for 5 min; 40 cycles of 95 °C for 30 s and 60 °C for 30 s. For each treatment, standard curves were produced based on a fivefold dilution series of cDNA as a template according to the linear regression model. The fixed threshold in this study was set to 500 to obtain all the threshold cycle (Ct) values of tested candidate reference genes. The qRT-PCR analyses were completed with three biological replicates and three technical replicates.Table 3 Functions, primer sequences, and amplicon characteristics of the candidate reference genes analyzed in this study. Gene symbol Gene name Gene ID (Putative) Function Primer sequences(5′-3′) aL (bp) bE (%) cR2 Actin Actin TR9961|c1_g1 Cytoskeletal structural protein F: CCATGTACCCTGGTATTGC 234 1.106 0.9984 R: TGTGGGAGGTGATGACTTA GAPDH glyceraldehyde-3-phosphate dehydrogenase TR3352|c0_g1 Glycolytic enzyme F: GGATTACCGACGCTACGC 232 0.977 0.9839 R: CGCACGCACAAGGATTTA NADH NADH dehydrogenase TR12676|c0_g1 Enzyme involved in redox reactions F: GTCAAACCTGGTGGCTAAA 182 0.941 0.9973 R: AGTCGTGGCGTCCATACAG AK arginine kinase TR3122|c0_g1 Key enzyme for cellular energy metabolism F: AGTACATAATTTCTACGAGGGT 169 1.014 0.9824 R: GACATGCCAGTTAAGGGA SDHB succinate dehydrogenase B TR11034|c0_g1 protein subunits of succinate dehydrogenase F: TCACGCCAGATTACCG 221 0.888 0.9998 R: TAGCTCCATGAACAGAAG RPL8 ribosomal protein L18 TR12462|c0_g1 Structural constituent of ribosome F: CCACAACCCAGACTCCA 179 0.935 0.9998 R: TAGGCCAGCAATTACGC RPL4 ribosomal protein L4 TR996|c0_g1 Structural constituent of ribosome F: AAAGCACCCATCAGACC 155 0.928 0.9961 R: CGGACACGAGGAATACG HSP68 heat shock protein68 TR7632|c0_g3 Molecular chaperone F: AAACGGGCTCGGGACA 245 0.955 0.9983 R: TCGACGGCGGGTGATA 18S 18S ribosomal RNA KT204362.1 Structural constituent of ribosome F: CGATGATGACGACGTGGTAGT 411 0.904 0.999 R: ACTACCACGTCGTCATCATCG EF1A elongation factor 1 a DQ005156.1 Catalysation of GTP-dependent binding of amynoayl-total RNA to the ribosome F: GGAACACGCTCTATTGGC 526 0.924 0.9989 R: CACGACCTACTGGGACTG aAmplicon length. bqRT-PCR efficiency (based on a standard curve). cReproducibility of the qRT-PCR. Data analysis The stability of the 10 candidate reference housekeeping genes was evaluated with the geNorm40, NormFinder73, and BestKeeper74 algorithms and the comparative delta Ct method75. Finally, we compared and ranked the tested candidate reference genes with the web-based RefFinder analytical tool (https://www.heartcure.com.au/for-researchers). Supplementary information Supplementary Table. Publisher's note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. These authors contributed equally: Xinan Li and Peipan Gong. Supplementary Information The online version contains supplementary material available at 10.1038/s41598-020-78974-z. Acknowledgements This study was supported by the National Key Research and Development Program of China (2018YFD0200501, 2017YFD0201703, 2016YFD0300705) and China Agriculture Research System (Award Number: CARS-3). Author contributions X.L., P.G., J.J. and X.Z. conceived and designed the research. X.L., P.G., M.L. and B.W. conducted the experiments. X.L., H.G., C.W. and X.Z. analyzed the data. X.L. and X.Z. wrote the manuscript. Y.Z., X.L. and L.W. revised the manuscript. All authors have read and approved the manuscript. Competing interests The authors declare no competing interests. ==== Refs References 1. Nolan T Hands RE Bustin SA Quantification of mRNA using real-time RT-PCR Nat. Protoc. 2006 1 1559 1582 10.1038/nprot.2006.236 17406449 2. Strube C Buschbaum S Wolken S Schnieder T Evaluation of reference genes for quantitative realtime PCR to investigate protein disulfide isomerase transcription pattern in the bovine lungworm Dictyocaulus viviparus Gene 2008 425 36 43 10.1016/j.gene.2008.08.001 18761062 3. Bustin SA Benes V Nolan T Pfaffl MW Quantitative real-time RT-PCR–a perspective J. Mol. Endocrinol. 2005 34 597 601 10.1677/jme.1.01755 15956331 4. Stephan L Tilmes V Hulskamp M Selection and validation of reference genes for quantitative Real-Time PCR in Arabis alpina PLoS ONE 2019 14 e0211172 10.1371/journal.pone.0211172 30830921 5. Chao WS Wang H Horvath DP Anderson JV Selection of endogenous reference genes for qRT-PCR analysis in Camelina sativa and identification of FLOWERING LOCUS C allele-specific markers to differentiate summer- and winter-biotypes Ind. Crops Prod. 2019 129 495 502 10.1016/j.indcrop.2018.12.017 6. Vandesompele J Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes Genome Biol. 2002 10.1186/gb-2002-3-7-research0034 12184808 7. Van Hiel MB Identification and validation of housekeeping genes in brains of the desert locust Schistocerca gregaria under different developmental conditions BMC Mol. Biol. 2009 10 56 10.1186/1471-2199-10-56 19508726 8. Pitino M Coleman AD Maffei ME Ridout CJ Hogenhout H Silencing of aphid genes by dsRNA feeding from plants PLoS ONE 2011 6 e25709 10.1371/journal.pone.0025709 21998682 9. Silva AX Jander G Samaniego H Ramsey JS Figueroa CC Insecticide resistance mechanisms in the green peach aphid Myzus persicae (Hemiptera: Aphididae) I: A transcriptomic survey PLoS ONE 2012 7 e36366 10.1371/journal.pone.0036366 22685538 10. Fu W Exploring valid reference genes for quantitative real-time PCR analysis in Plutella xylostella (Lepidoptera: Plutellidae) Int. J. Biol. Sci. 2013 9 792 802 10.7150/ijbs.5862 23983612 11. Yang C Pan H Liu Y Zhou X Selection of reference genes for expression analysis using quantitative real-time PCR in the pea aphid, Acyrthosiphon pisum (Harris) (Hemiptera, Aphidiae) PLoS ONE 2014 9 e110454 10.1371/journal.pone.0110454 25423476 12. Radonić A Guideline to reference gene selection for quantitative real-time PCR Biochem. Biophys. Res. Commun. 2004 313 856 862 10.1016/j.bbrc.2003.11.177 14706621 13. Dheda K The implications of using an inappropriate reference gene for real-time reverse transcription PCR data normalization Anal. Biochem. 2005 344 141 143 10.1016/j.ab.2005.05.022 16054107 14. Dong J Li J Huang J Wu J Identification of suitable reference genes for miRNA quantitation in bumblebee (Hymenoptera: Apidae) response to reproduction Apidologie 2019 50 40 50 10.1007/s13592-018-0616-9 15. Yang X Pan H Yuan L Zhou X Reference gene selection for RT-qPCR analysis in Harmonia axyridis , a global invasive lady beetle Sci. Rep. 2018 8 2689 10.1038/s41598-018-20612-w 29426915 16. Qu C Selection and evaluation of reference genes for expression analysis using quantitative real-time PCR in the Asian Ladybird Harmonia axyridis (Coleoptera: Coccinellidae) PLoS ONE 2018 13 e0192521 10.1371/journal.pone.0192521 29889877 17. Lü J Selection of appropriate reference genes for RT-qPCR analysis in Propylea japonica (Coleoptera: Coccinellidae) PLoS ONE 2018 13 e0208027 10.1371/journal.pone.0208027 30481225 18. Yang C Pan H Liu Y Zhou X Temperature and development impacts on housekeeping gene expression in cowpea aphid, Aphis craccivora (Hemiptera: Aphidiae) PLoS ONE 2015 10 e0130593 10.1371/journal.pone.0130593 26090683 19. Lü J Selection and validation of reference genes for RT-qPCR analysis of the ladybird beetle Henosepilachna vigintioctomaculata Front. Physiol. 2018 9 1614 10.3389/fphys.2018.01614 30483159 20. Xu J Lu M Cui Y Du Y Selection and evaluation of reference genes for expression analysis using qRT-PCR in Chilo suppressalis (Lepidoptera: Pyralidae) J. Econ. Entomol. 2017 110 683 691 10.1093/jee/tow297 28115499 21. Tan Y Zhou XR Pang BP Reference gene selection and evaluation for expression analysis using qRT-PCR in Galeruca daurica (Joannis) Bull. Entomol. Res. 2017 107 359 368 10.1017/s0007485316000948 27819206 22. Chang Y Selection and validation of reference genes for quantitative real-time PCR analysis under different experimental conditions in the leafminer Liriomyza trifolii (Diptera: Agromyzidae) PLoS ONE 2017 12 e0181862 10.1371/journal.pone.0181862 28746411 23. Arya SK Reference genes validation in Phenacoccus solenopsis under various biotic and abiotic stress conditions Sci. Rep. 2017 7 13520 10.1038/s41598-017-13925-9 29051594 24. Yang C Selection of reference genes for RT-qPCR analysis in Coccinella septempunctata to assess un-intended effects of RNAi transgenic plants Front. Plant Sci. 2016 7 1672 10.3389/fpls.2016.01672 27877186 25. Koramutla MK Aminedi R Bhattacharya R Comprehensive evaluation of candidate reference genes for qRT-PCR studies of gene expression in mustard aphid, Lipaphis erysimi (Kalt) Sci. Rep. 2016 6 25883 10.1038/srep25883 27165720 26. Kang Z Evaluation of the reference genes for expression analysis using quantitative real-time polymerase chain reaction in the green peach aphid, Myzus persicae Insect Sci. 2017 24 222 234 10.1111/1744-7917.12310 26749166 27. Cristiano G Validation of reference genes for qRT-PCR analysis in Megoura viciae (Hemiptera Aphididae) Bull. Insectol. 2016 69 229 238 28. Rabbinge R Mantel WP Monitoring for cereal aphids in winter wheat Neth. J. Plant Pathol. 1981 87 25 29 10.1007/bf01981397 29. Cannon RJC Summer populations of the cereal aphid Metopolophium dirhodum (Walker) on winter wheat: Three contrasting years J. Appl. Ecol. 1986 23 101 114 10.2307/2403084 30. Honěk A Factors determining the peak abundance of Metopolophium dirhodum (Homoptera: Aphididae) on cereals Bull. Entomol. Res. 1991 81 57 64 10.1017/s0007485300053244 31. Ma C Hau B Poehling H Effects of pattern and timing of high temperature exposure on reproduction of the rose grain aphid, Metopolophium dirhodum Entomol. Exp. Appl. 2004 110 65 71 10.1111/j.0013-8703.2004.00123.x 32. Zuñiga E Control biológico de los afidos de los cereales en Chile. I. Revisión histórica y líneas de trabajo Agric. Tec. 1986 46 475 477 33. Honek A Martinkova Z Saska P Dixon AFG Aphids (Homoptera: Aphididae) on winter wheat: Predicting maximum abundance of Metopolophium dirhodum J. Econ. Entomol. 2018 10.1093/jee/toy157 29897464 34. Holt J Griffiths E Wratten SD The influence of wheat growth stage on yield reductions caused by the rose-grain aphid, Metopolophium dirhodum Ann. Appl. Biol. 1984 105 7 14 10.1111/j.1744-7348.1984.tb02797.x 35. Kennedy TF Connery J Grain yield reductions in spring barley due to barley yellow dwarf virus and aphid feeding Ir. J. Agric. Food Res. 2005 44 111 128 36. Ross DT Systematic variation in gene expression patterns in human cancer cell lines Nat. Genet. 2000 24 227 235 10.1038/73432 10700174 37. Solanas M Moral R Escrich E Unsuitability of using ribosomal RNA as loading control for Northern blot analyses related to the imbalance between messenger and ribosomal RNA content in rat mammary tumors Anal. Biochem. 2001 288 99 102 10.1006/abio.2000.4889 11141312 38. Bustin SA The MIQE guidelines: Minimum information for publication of quantitative real-time PCR experiments Clin. Chem. 2009 55 611 622 10.1373/clinchem.2008.112797 19246619 39. Mao J Zeng F Plant-mediated RNAi of a gap gene-enhanced tobacco tolerance against the Myzus persicae Transgenic Res. 2014 23 145 152 10.1007/s11248-013-9739-y 23949691 40. Vandesompele J Paepe AD Speleman F Elimination of primer-dimer artifacts and genomic coamplification using a two-step SYBR green I real-time RT-PCR Anal. Biochem. 2002 303 95 98 10.1006/abio.2001.5564 11906156 41. Bustin SA Quantification of mRNA using real-time reverse transcription PCR (RT-PCR): Trends and problems J. Mol. Endocrinol. 2002 29 23 29 10.1677/jme.0.0290023 12200227 42. Gutierrez L The lack of a systematic validation of reference genes: A serious pitfall undervalued in reverse transcription-polymerase chain reaction (RT-PCR) analysis in plants Plant Biotechnol. J. 2008 6 609 618 10.1111/j.1467-7652.2008.00346.x 18433420 43. Liang P Guo Y Zhou X Gao X Expression profiling in Bemisia tabaci under insecticide treatment: Indicating the necessity for custom reference gene selection PLoS ONE 2014 9 e87514 10.1371/journal.pone.0087514 24498122 44. Yang C Selection of reference genes for RT-qPCR analysis in a predatory biological control agent, Coleomegilla maculata (Coleoptera: Coccinellidae) Sci. Rep. 2015 5 18201 10.1038/srep18201 26656102 45. Pan H Yang X Siegfried BD Zhou X A comprehensive selection of reference genes for RT-qPCR analysis in a predatory lady beetle, Hippodamia convergens (Coleoptera: Coccinellidae) PLoS ONE 2015 10 e0125868 10.1371/journal.pone.0125868 25915640 46. Li R Reference gene selection for qRT-PCR analysis in the sweetpotato whitefly, Bemisia tabaci (Hemiptera: Aleyrodidae) PLoS ONE 2013 8 e53006 10.1371/journal.pone.0053006 23308130 47. Zhu X Selection and evaluation of reference genes for expression analysis using qRT-PCR in the beet armyworm Spodoptera exigua (Hübner)(Lepidoptera: Noctuidae) PLoS ONE 2014 9 e84730 10.1371/journal.pone.0084730 24454743 48. Thellin O Housekeeping genes as internal standards: Use and limits J. Biotechnol. 1999 75 291 295 10.1016/S0168-1656(99)00163-7 10617337 49. Lü J Yang C Zhang Y Pan H Selection of reference genes for the normalization of RT-qPCR data in gene expression studies in insects: A systematic review Front. Physiol. 2018 9 1560 10.3389/fphys.2018.01560 30459641 50. Shi C Evaluation of housekeeping genes for quantitative real-time PCR analysis of Bradysia odoriphaga (Diptera: Sciaridae) Int. J. Mol. Sci. 2016 17 1034 10.3390/ijms17071034 51. Shakeel M Zhu X Kang T Wan H Li J Selection and evaluation of reference genes for quantitative gene expression studies in cotton bollworm, Helicoverpa armigera (Lepidoptera: Noctuidae) J. Asia-Pac. Entomol. 2015 18 123 130 10.1016/j.aspen.2015.01.001 52. De Boer ME Reference genes for QRT-PCR tested under various stress conditions in Folsomia candida and Orchesella cincta (Insecta, Collembola) BMC Mol. Biol. 2009 10 54 10.1186/1471-2199-10-54 19486513 53. Frost P Nilsen F Validation of reference genes for transcription profiling in the salmon louse, Lepeophtheirus salmonis , by quantitative real-time PCR Vet. Parasitol. 2003 118 169 174 10.1016/j.vetpar.2003.09.020 14651887 54. Sun W Jin Y He L Lu WC Li M Suitable reference gene selection for different strains and developmental stages of the carmine spider mite, Tetranychus cinnabarinus , using quantitative real-time PCR J. Insect Sci. 2010 10 208 10.1673/031.010.20801 21265619 55. Hornakova D Matouskova P Kindl JÍ Valterova I Pichova I Selection of reference genes for real-time polymerase chain reaction analysis in tissues from Bombus terrestris and Bombus lucorum of different ages Anal. Biochem. 2010 397 118 120 10.1016/j.ab.2009.09.019 19751695 56. An X Hou M Liu Y Reference gene selection and evaluation for gene expression studies using qRT-PCR in the white-backed planthopper, Sogatella furcifera (Hemiptera: Delphacidae) J. Econ. Entomol. 2015 109 879 886 10.1093/jee/tov333 57. Chapuis MP Assessment and validation of a suite of reverse transcription-quantitative PCR reference genes for analyses of density-dependent behavioural plasticity in the Australian plague locust BMC Mol. Biol. 2011 12 7 10.1186/1471-2199-12-7 21324174 58. Sun M Lu MX Tang XT Du YZ Exploring valid reference genes for quantitative real-time PCR analysis in Sesamia inferens (Lepidoptera: Noctuidae) PLoS ONE 2015 10 e0115979 10.1371/journal.pone.0115979 25585250 59. Zhai Y Identification and validation of reference genes for quantitative real-time PCR in Drosophila suzukii (Diptera: Drosophilidae) PLoS ONE 2014 9 e106800 10.1371/journal.pone.0106800 25198611 60. Mamidala P Rajarapu SP Jones SC Mittapalli O Identification and validation of reference genes for quantitative real-time polymerase chain reaction in Cimex lectularius J. Med. Entomol. 2011 48 947 951 10.1603/me10262 21845960 61. Bustin SA Absolute quantification of mRNA using real-time reverse transcription polymerase chain reaction assays J. Mol. Endocrinol. 2000 25 169 193 10.1677/jme.0.0250169 11013345 62. Bagnall NH Kotze AC Evaluation of reference genes for real-time PCR quantification of gene expression in the Australian sheep blowfly, Lucilia cuprina Med. Vet. Entomol. 2010 24 176 181 10.1111/j.1365-2915.2010.00866.x 20604863 63. Majerowicz D Looking for reference genes for real-time quantitative PCR experiments in Rhodnius prolixus (Hemiptera: Reduviidae) Insect Mol. Biol. 2011 20 713 722 10.1111/j.1365-2583.2011.01101.x 21929722 64. Paim RM Validation of reference genes for expression analysis in the salivary gland and the intestine of Rhodnius prolixus (Hemiptera, Reduviidae) under different experimental conditions by quantitative real-time PCR BMC Res. Notes 2012 5 128 10.1186/1756-0500-5-128 22395020 65. Maroniche GA Sagadín M Mongelli VC Truol GA del Vas M Reference gene selection for gene expression studies using RT-qPCR in virus-infected planthoppers Virol. J. 2011 8 308 10.1186/1743-422x-8-308 21679431 66. Shen G Jiang H Wang X Wang J Evaluation of endogenous references for gene expression profiling in different tissues of the oriental fruit fly Bactrocera dorsalis (Diptera: Tephritidae) BMC Mol. Biol. 2010 11 76 10.1186/1471-2199-11-76 20923571 67. Tricarico C Quantitative real-time reverse transcription polymerase chain reaction: Normalization to rRNA or single housekeeping genes is inappropriate for human tissue biopsies Anal. Biochem. 2002 309 293 300 10.1016/s0003-2697(02)00311-1 12413463 68. Zhao L Jones WA Expression of heat shock protein genes in insect stress responses Invertebr. Surv. J. 2012 9 93 101 10.1155/2012/484919 69. Veazey KJ Golding MC Selection of stable reference genes for quantitative rt-PCR comparisons of mouse embryonic and extra-embryonic stem cells PLoS ONE 2011 6 e27592 10.1371/journal.pone.0027592 22102912 70. Ling D Salvaterra PM Robust RT-qPCR data normalization: Validation and selection of internal reference genes during post-experimental data analysis PLoS ONE 2011 6 e17762 10.1371/journal.pone.0017762 21423626 71. Zuo Y Regional susceptibilities of Rhopalosiphum padi (Hemiptera: Aphididae) to ten insecticides Fla. Entomol. 2016 99 269 275 10.1653/024.099.0217 72. Wilkinson TL Ishikawa H On the functional significance of symbioticm icroorganisms in the Homoptera: A comparative study of Acyrthosiphon pisum and Nilaparvata lugens Physiol. Entomol. 2001 26 86 93 73. Andersen CL Jensen JL Ørntoft TF Normalization of real-time quantitative reverse transcription-PCR data: A model-based variance estimation approach to identify genes suited for normalization, applied to bladder and colon cancer data sets Can. Res. 2004 64 5245 5250 10.1158/0008-5472.CAN-04-0496 74. Pfaffl MW Tichopad A Prgomet C Neuvians TP Determination of stable housekeeping genes, differentially regulated target genes and sample integrity: BestKeeper—Excel-based tool using pair-wise correlations Biotechnol. Lett. 2004 26 509 515 10.1023/B:BILE.0000019559.84305.47 15127793 75. Silver N Best S Jiang J Thein SL Selection of housekeeping genes for gene expression studies in human reticulocytes using real-time PCR BMC Mol. Biol. 2006 7 33 10.1186/1471-2199-7-33 17026756