Identification and Validation of Reference Genes for Expression Analysis Using qRT-PCR in Cimex hemipterus (Hemiptera: Cimicidae)

Simple Summary Reference genes are the basis for exploring the function of target genes in insects. Recently, we identified a large number of chemoreceptor genes by transcriptome sequencing in different tissues of Cimex hemipterus (F.) (Hemiptera: Cimicidae). Additionally, we plan to study the expression patterns of these genes and clarify their functions. However, there is no report on the screening and application of the reference genes in C. hemipterus. In the present study, the expression-level stability of 10 candidate reference genes in C. hemipterus under various experimental conditions was evaluated through qRT-PCR combined with RefFinder (which integrates four computational programs: geNorm, NormFinder, BestKeeper, and ∆Ct). Results show that the optimal combination of reference genes for each experimental condition was as follows: RPL8 and EF1α for the developmental stage, RPL8 and RPS16 for adult sex, RPL8 and RPL11 for adult tissue, RPL8 and β-tubulin for gas stimulation, and RPL8 and NADH for temperature. The results will lay a foundation for the follow-up study on the expression and function of the target genes of tropical bed bugs. Abstract The tropical bed bug, Cimex hemipterus (F.) (Hemiptera: Cimicidae) is an important public-health pest that feeds on the blood of humans and some other animals. To explore the function of the target genes of C. hemipterus, it is essential to select suitable reference genes for the accurate quantification of gene expression. Here, we selected 10 frequently used reference genes in insects and evaluated their stability in C. hemipterus under various biotic (developmental stage, sex, and tissue) and abiotic (gas stimulation and temperature) conditions through RefFinder (which integrates four computational programs: geNorm, NormFinder, BestKeeper, and ∆Ct). Results indicate that the optimal combination of reference genes for each experimental condition was as follows: RPL8 and EF1α for the developmental stage (eggs, early instar nymphs, late instar nymphs, and adults), RPL8 and RPS16 for adult sex, RPL8 and RPL11 for adult tissue (head, thorax, abdomen, and legs), RPL8 and β-tubulin for gas stimulation (air and carbon dioxide), and RPL8 and NADH for temperature (0, 5, 17, 30, and 37 °C). Finally, the expression pattern of the HSP70 and GR21 genes were analyzed, and the results highlight the importance of appropriate reference-gene selection. Our results provide a comprehensive list of optimal reference genes from C. hemipterus for the first time, which will contribute to accurately analyzing the expression of target genes.

The Cimex hemipterus used in this study were collected from Guangzhou City, Guangdong Province, China in 2018 [19]. They were maintained in artificial climate chambers (26 ± 1 • C, 40 ± 10% relative humidity (RH) with a photoperiod of 12:12 (L:D) h) according to the methods of   [35].

Experimental Condition and Sample Collection
A total of five experimental conditions (developmental stage, sex, tissue, gas stimulation, and temperature) were set up to assess the stability of the 10 candidate reference genes in C. hemipterus. For each replication, the developmental-stage sample includes 50 eggs, 20 1st-2nd instar nymphs, 8 4th-5th instar nymphs, and 5 male adults; the sex sample contains 5 male adults and 5 female adults; the body tissue sample includes the head (containing antennae and mouthparts), thorax, abdomen, and legs, which were dissected from 50 male adults; the gas stimulation sample contains 5 CO 2 -treated (the concentration of CO 2 was 1.6%) male adults and 5 air-treated male adults-all the male adults were stimulated by CO 2 or air for 3 h with a flow rate of 500 mL/min, respectively; the temperature treatment sample contains male adults that were exposed to 0, 5, 17, 30, and 37 • C for 2 h, with 4 individuals per temperature treatment. Each experimental condition was replicated three times. The collected samples were immediately frozen with liquid nitrogen and temporarily stored under the condition of −80 • C until RNA extraction.

RNA Extraction and cDNA Synthesis
TRIzol reagent (Takara Biomedical Technology (Beijing) Co., Ltd., Beijing, China) was used for extracting total RNA from each sample following the manufacturer's manual. The concentration and purity of RNAs were determined by a NanoDrop One Spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA). The OD260/280 value of the total RNA was between 1.9 and 2.1, which was considered to meet the quality requirements for subsequent experiments. Then, the integrity of RNAs was examined by 0.8% agarose gel electrophoresis. First-strand cDNAs were synthesized by reverse transcribing 1 µg of total RNA with the PrimeScript RT Kit (Takara Biomedical Technology (Beijing) Co., Ltd., Beijing, China).
The 20 µL reaction system included 10 µL 2 × SYBR Green Premix (TIANGEN, China), 0.6 µL of 10 µM primers, 6.8 µL ddH 2 O, and 2 µL diluted cDNA. Reactions were conducted using a CFX-96 real-time PCR system (Bio-Rad, Hercules, CA, USA). The amplification program was as follows: 95 • C for 10 min, 40 cycles of 95 • C for 10 s, and 60 • C for 32 s. The presence of a single peak in the melting curve analysis was selected to confirm the gene-specific amplification, as well as to get rid of the non-specific amplification and the primer-dimer generation.
The amplification efficiency and correlation coefficient of each candidate reference gene were measured by applying the slope of the standard curve generated by the serial five-fold diluted cDNA template [1]. The amplification efficiency (E) of each primer pair was calculated according to the equation [5]:

Determination of Reference Gene Expression Stability
The expression stability of the 10 candidate reference genes was assessed using RefFinder (https://www.heartcure.com.au/reffinder/, accessed on 10 July 2022), which integrates four major computational programs, i.e., geNorm [36], NormFinder [7], Best-Keeper [4], and the ∆Ct method [37]. The suitable number of reference genes for normalizing target gene expression was determined by computing the pairwise variation (V) using geNorm [36]. If the V-value (Vn/Vn+1) is lower than 0.15, it indicates using n genes is enough for target gene expression normalization.

Validation of the Candidate Reference Genes
To further verify the accuracy of our experimental results, the expression profiles of the heat shock protein (HSP70) gene and gustatory receptor 21 (GR21) gene in C. hemipterus were normalized by the top two most stable (RPL8 and RPL11) and least stable (RPL13 and RPS16) candidates in the head (including antennae and mouthparts), thorax, abdomen, and legs. The forward primer of HSP70 was 5 -GAACGGGAACCAAGGAGGC-3 , and the reverse primer of HSP70 was 5 -TCGACGGCGAATCTCAGCA-3 . The forward primer of GR21 was 5 -CGGGCATGGTGTTTATCAAG-3 , and the reverse primer of GR21 was 5 -CTCCCAGTGTAACTCGAGGG-3 . The relative expression levels of HSP70 and GR21 in different tissues (considered as different treatments) were computed using the 2 −∆∆Ct method [38]. Each treatment was replicated 3-5 times.
One-way analysis of variance (ANOVA) was used to analyze the expression differences between HSP70 and GR21 in different tissues, and Tukey's honest significant difference (HSD) test was selected for separating the means. Statistical analysis was conducted using SPSS version 22.0 [39].

Primer Specificity and Efficiency
Melting curves of all primer pairs were a single peak, and no nonspecific amplification products were observed, indicating that the primers were highly specific. The amplification efficiency (E) of the primers for the 10 reference genes ranged from 87.11% for NADH to 111.07% for RPL8, and the correlation coefficients varied from 0.9806 to 0.9990 (Table 1).

Ct Values of Candidate Reference Genes
The cycle threshold (Ct) values of the 10 reference genes under five experimental conditions were used for evaluating their expression levels. As shown in Figure 1 and One-way analysis of variance (ANOVA) was used to analyze the expression differences between HSP70 and GR21 in different tissues, and Tukey's honest significant difference (HSD) test was selected for separating the means. Statistical analysis was conducted using SPSS version 22.0 [39].

Primer Specificity and Efficiency
Melting curves of all primer pairs were a single peak, and no nonspecific amplification products were observed, indicating that the primers were highly specific. The amplification efficiency (E) of the primers for the 10 reference genes ranged from 87.11% for NADH to 111.07% for RPL8, and the correlation coefficients varied from 0.9806 to 0.9990 (Table 1).

Ct Values of Candidate Reference Genes
The cycle threshold (Ct) values of the 10 reference genes under five experimental conditions were used for evaluating their expression levels. As shown in Figure 1 and Table 2, the Ct values ranged from 14.92 to 32.34, and the most values were distributed between 18 and 27. Actin was the most abundant transcript (the mean Ct value was 16.46 ± 0.40), and NADH was the least abundant candidate reference gene (the mean Ct value was 24.61 ± 0.56) ( Table 2).

Stability of the 10 Reference Genes under Different Experimental Conditions
For the developmental stage, there was a significant difference between the two most stable reference genes recommended by the four algorithms (Table 3). In the geNorm analysis, RPL8 and RPL11 had the best stability. In the NormFinder analysis, EF1α and NADH were the two most stable genes. The BestKeeper analysis showed that NADH and β-tubulin were the top two ranked genes. Additionally, in the ∆Ct analysis, RPL8 and EF1α were identified as the two most stable genes. Based on the RefFinder analysis results, RPL8 was the most stable reference gene in the developmental stage condition, followed by EF1α, NADH, RPL11, GAPDH, RPS16, β-tubulin, RPL13, α-tubulin, and Actin ( Figure 2). The pairwise variation value of V2/3 calculated by geNorm was below the threshold value of 0.15 ( Figure 3). Thus, in different developmental stages, RPL8 and EF1α were the suitable combination of reference genes.
Insects 2022, 13, x FOR PEER REVIEW 7 of 15 followed by β-tubulin, Actin, RPL11, EF1α, α-tubulin, GAPDH, NADH, RPL13, and RPS16 ( Figure 2). The geNorm analysis results indicated that the pair-wise value of V2/3 was lower than 0.15 (Figure 3), suggesting that the combination of RPL8 and β-tubulin was sufficient for reliable normalization in different gas stimulation samples. For temperature treatment, the top two most stable genes from the four computational programs differed from each other. They were RPL8 and RPL13 by geNorm, RPL8, and β-tubulin by NormFinder, NADH and Actin by BestKeeper, and RPL8 and NADH by ΔCt (Table 3). Based on the RefFinder analysis, the most stable reference gene under the conditions of different temperatures was RPL8, followed by NADH, βtubulin, RPL13, EF1α, GAPDH, RPL11, Actin, α-tubulin, and RPS16 ( Figure 2). The pairwise variation value of V2/3 calculated by the geNorm analysis program was less than 0.15 (Figure 3), indicating that the combination of RPL8 and NADH was adequate for optimal normalization in different temperature treatment samples.     For bed bug sex, the two most stable genes from the four computational programs differ significantly. For example, RPS16 and EF1α were determined to be the top two most stably expressed referenced genes by geNorm, RPL8 and RPL13 for NormFinder, RPS16 and RPL11 for BestKeeper, and RPL8 and RPL11 for ∆Ct (Table 3). Based on the integrated analysis results of RefFinder, the most stable reference gene in different sex conditions was RPL8, followed by RPS16, RPL11, EF1α, GAPDH, RPL13, α-tubulin, Actin, β-tubulin, and NADH (Figure 2). The results of the geNorm analysis indicate that the pairwise variation value of V2/3 was less than the threshold value of 0.15 ( Figure 3). Therefore, the combination of RPL8 and RPS16 was recommended for normalizing the qRT-PCR data in the different sexes of adult samples by geNorm.
For bed bug tissue, geNorm and NormFinder had the same top two genes in stability ranking: RPL8 and RPL11 (Table 3). However, the top two in BestKeeper were RPL8 and α-tubulin, and the top two in ∆Ct were EF1α and RPL8 ( Table 3). The RefFinder analysis showed that RPL8 was the most stable reference gene under this experimental condition, followed by RPL11, EF1α, α-tubulin, GAPDH, Actin, β-tubulin, NADH, RPL13, and RPS16 ( Figure 2). The pairwise variation value of V2/3 was less than the cut-off value of 0.15 according to the geNorm analysis. (Figure 3). Therefore, RPL8 and RPL11 were recognized as suitable reference genes for normalizing qRT-PCR data in the different tissue samples.
For gas stimulation, NormFinder and ∆Ct had the same top two gene stability rankings: RPL8 and β-tubulin (Table 3). However, based on geNorm, Actin and β-tubulin were the top two most stable genes. EF1α and GAPDH were determined to be the top two most stably expressed genes by BestKeeper (Table 3). Based on the analysis results of RefFinder, the most stable gene under the treatment of gas stimulation was RPL8, followed by βtubulin, Actin, RPL11, EF1α, α-tubulin, GAPDH, NADH, RPL13, and RPS16 ( Figure 2). The geNorm analysis results indicated that the pair-wise value of V2/3 was lower than 0.15 (Figure 3), suggesting that the combination of RPL8 and β-tubulin was sufficient for reliable normalization in different gas stimulation samples.
For temperature treatment, the top two most stable genes from the four computational programs differed from each other. They were RPL8 and RPL13 by geNorm, RPL8, and β-tubulin by NormFinder, NADH and Actin by BestKeeper, and RPL8 and NADH by ∆Ct (Table 3). Based on the RefFinder analysis, the most stable reference gene under the conditions of different temperatures was RPL8, followed by NADH, β-tubulin, RPL13, EF1α, GAPDH, RPL11, Actin, α-tubulin, and RPS16 ( Figure 2). The pairwise variation value of V2/3 calculated by the geNorm analysis program was less than 0.15 (Figure 3), indicating that the combination of RPL8 and NADH was adequate for optimal normalization in different temperature treatment samples.

Validation of Candidate Reference Genes in Different Tissues
The expression trend of HSP70 in different tissues was similar in the two combinations, with higher expression in the thorax and abdomen, and lower expression in the head and legs. When the two most stable reference genes (RPL8 and RPL11) were selected as reference genes, the relative expression levels of HSP70 in the legs were significantly higher than in the head (F 3,16 = 17.1, p < 0.0001; Figure 4). However, when the two least stable reference genes (RPL13 and RPS16) were selected, no significant difference was observed in the relative expression of HSP70 between the head and legs (F 3,16 = 9.1, p = 0.001; Figure 4). For another gene, the relative expression trend of GR21 in different tissues was inconsistent in the two reference-gene combinations (RPL8 and RPL11, RPL13 and RPS16). The relative expression of GR21 in the head was the highest, and the relative expression in the four tissues was significantly different, when normalized by RPL8 and RPL11 (F 3,8 = 439.9, p < 0.0001; Figure 5). While, when normalized by RPL13 and RPS16, the tissue with the highest relative expression of GR21 was the thorax. However, there was no significant difference in the relative expression among the head, thorax, and legs; in addition, no significant difference was observed among the head, abdomen, and legs ( Figure 5). Therefore, these results further indicate that the expression levels of the target genes will be different when normalized by different reference-gene combinations.

Discussion
In the present study, the expression-level stability of 10 candidate reference genes in C. hemipterus under various experimental conditions was evaluated through qRT-PCR. According to the final evaluation results of RefFinder, it was found that the stable expressed reference genes recommended for different developmental stages, different tissues of the adult males, different sexes, different temperature treatments, and bed bugs stimulated by different gases varied. Among them, RPL8 had the least variable expression levels under the five experimental conditions, which was consistent with the report of Zhu Figure 5. Tissue expression patterns of GR21 gene in C. hemipterus. The relative expression levels of GR21 in the head (including antennae and mouthparts), thorax, abdomen, and legs were normalized to the top two most stable (A, RPL8 and RPL11) and least stable (B, RPL13 and RPS16) candidate reference genes, respectively. Different lower letters above bars represent significant differences in gene expression among different tissues (p < 0.05, Tukey's HSD test).

Discussion
In the present study, the expression-level stability of 10 candidate reference genes in C. hemipterus under various experimental conditions was evaluated through qRT-PCR. According to the final evaluation results of RefFinder, it was found that the stable expressed reference genes recommended for different developmental stages, different tissues of the adult males, different sexes, different temperature treatments, and bed bugs stimulated by different gases varied. Among them, RPL8 had the least variable expression levels under the five experimental conditions, which was consistent with the report of Zhu et al. [34]. In the research of Zhu et al. [34], the stability of four reference genes (RPL11, RPL8, RPS16, and HSP70) under the conditions of different populations, different developmental stages, different tissues, and dsRNA injected in C. lectularius was compared; that study also identified RPL8 as an optimal reference gene. We infer that this consistency may be related to the close genetic relationship between C. hemipterus and C. lectularius. Another report about the selection of reference genes in C. lectularius showed that RPL18 was defined as the best reference gene for various tissues and developmental stages of pesticide-exposed and susceptible strains, and their candidate reference genes include Actin, β-tubulin, EF1α, GAPDH, syntaxin (SYN), RPL18, and ubiquitin-C (UBC) [33]. However, our study did not include RPL18 as a candidate reference gene. Whether there was a difference in the stability of RPL8 and RPL18 in C. hemipterus requires further study.
Actin and GAPDH, the most widely used internal reference genes in animals, insects, and plants, did not show good expression stability in C. hemipterus (Figure 2, Table 3). This significant variance in the expression stability has been reported in many previous studies [6,8,9,33,40,41]. EF1α [33,42] and NADH [43] also have been widely used as reference genes in different species and experimental conditions. In our experiments, EF1α was one of the top two most stable genes in the developmental-stage sample set, but it did not perform well among the other four factors. Similarly, NADH was one of the top two most stable genes in different temperature treatment samples; however, its stabilization was greatly reduced under the conditions of different developmental stages and sexes. These results further demonstrate that no single reference gene is suitable for all species and different treatment conditions.
The ribosomal protein (RP) is involved in translation and protein synthesis and is one of the most highly conserved proteins in all biological cells. The RP coding gene was considered to be one stable expressed reference gene and had been widely used in insect molecular research to regulate the level of gene expression over the past 10 years [44][45][46][47][48]. These genes include RPL8 under the conditions of different populations, different developmental stages, different tissues, dsRNA injected in C. lectularius [34], different geographic populations, across developmental stages, different wing dimorphisms, different temperature treatments, and different antibiotic treatments of Metopolophium dirhodum (Walker) (Hemiptera: Aphididae) [43]; RPL13 in the different tissues of Rhopalosiphum padi L. (Hemiptera: Aphididae) [49]; and RPL11 in the conditions of different developmental stages and temperatures in Aphis craccivora Koch (Hemiptera: Aphidiae) [50]. In this study, the RPL8 gene exhibited the most stable expression level in C. hemipterus under the five experimental conditions, RPL11 was recommended as one of the top two most stable expressed genes in different adult tissues, and RPS16 was also recommended for different sexes. Our results further expand the application scope of the RP coding gene as an internal reference gene in hemiptera insects.
It has been widely accepted that introducing two or more reference genes to normalize target-gene expression levels can reduce various errors and ensure the accuracy of results in experiments [51,52]. However, using a poor stable reference gene may yield the incorrect expression pattern of the target gene and lead to erroneous interpretations [50]. In the present experiment, the expression level of GR21 was inconsistent among different tissues when normalized to the top two most stable (RPL8 and RPL11) and least stable (RPL13 and RPS16) reference genes. Therefore, it is pivotal to choose and confirm the best stable reference genes for certain insect species under specific conditions. Over the past 20-30 yr period, there was a significant upward trend of bed bug infestations worldwide. The effective control of C. hemipterus has become a common issue in many infested countries or regions, such as China [19,20,53], southeast Asia [54,55], and South Africa [56]. Our results can provide the essential information for target-gene expression analysis and functional studies of this pest in order to lay a foundation for further exploration of new C. hemipterus control methods/strategies. However, this study only screened the reference genes of C. hemipterus under three biotic (developmental stage, tissue, and sex) and two abiotic (gas stimulation and temperature) conditions. Whether the results are also suitable for other experimental conditions (such as RNAi, insecticide treatment, etc.) in C. hemipterus or other species needs to be further verified.

Conclusions
In this study, we selected 10 frequently used reference genes in insects and evaluated their stability in C. hemipterus under various biotic (developmental stage, sex, and tissue) and abiotic (gas stimulation and temperature) conditions using RefFinder. The results show that the most suitable candidate combinations of reference genes were as follows: RPL8 and EF1α for the developmental stage (eggs, early instar nymphs, late instar nymphs, and adults), RPL8 and RPS16 for adult sex, RPL8 and RPL11 for adult tissue (head, thorax, abdomen, and legs), RPL8 and β-tubulin for gas stimulation (air and carbon dioxide), and RPL8 and NADH for temperature treatment. These findings will provide the basic information for subsequent target-gene-expression analysis and functional studies in C. hemipterus.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.