Selection and Evaluation of Reference Genes for qRT-PCR in Spodoptera frugiperda (Lepidoptera: Noctuidae)

Simple Summary The fall armyworm Spodoptera frugiperda (J. E. Smith) (Lepidoptera: Noctuidae) is an extremely important omnivorous agricultural pest, it poses a severe threat to food security and agricultural production. Quantitative real-time PCR (qRT-PCR) is an important molecular technology widely used for expression profile analyses of various target genes. It is essential to use reference genes as the benchmark to eliminate various errors and normalize the qRT-PCR analysis. In our study, 10 reference genes were evaluated under six experimental conditions, including developmental stages, tissues, mating status, hormones, diets, and temperatures. Finally, the expression profile of the target gene SfrOBP1 in various tissues of S. frugiperda was evaluated to verify the accuracy of the results. This study will provide a preliminary evaluation of reference genes of S. frugiperda, which can be beneficial to the further research of functional gene expression. Abstract As an accurate and convenient technique, the qRT-PCR is always used in the quantitative expression analysis of functional genes. Normalization of the data relies on stable reference genes. The fall armyworm Spodoptera frugiperda (J. E. Smith) is an important invasive and migratory pest that seriously threatens corn production around the world. In this paper, we selected 10 candidate reference genes (18S, AK, RPL10, RPS24, 28S, SOD, ATP, GAPDH, ACT, and a-TUB) and determined their expression levels under different conditions (different developmental stages, various tissues, mating status, hormones, diets, and temperatures). Subsequently, the stability of reference genes was evaluated by four algorithms (Delta Ct method, geNorm, NormFinder, BestKeeper). The optimal combination of reference genes for each treatment was obtained by geNorm. Finally, the comprehensive ranks were determined by the online tool RefFinder. Results showed that the most stable reference genes were SOD, RPL10, and RPS24 for developmental stages, α-TUB, RPL10, and ATP for different tissues, AK, RPL10, and 18S for mating status, 18S and AK under hormone treatment, 18S, RPL10, and SOD under diet treatment, RPL10, 18S, and RPS24 under temperature treatment. This study confirmed recent data on a few reference genes and provided an evaluation of a number of additional reference genes of S. frugiperda under various conditions.


Introduction
The fall armyworm Spodoptera frugiperda (J. E. Smith) (Lepidoptera: Noctuidae) is indigenous to the tropical and subtropical regions of the American continent, with high reproductive capacity, a wide range of hosts, and strong migration ability [1]. Therefore, this pest seriously threatens corn production and food security around the world [2][3][4]. At present, the S. frugiperda has invaded 46 countries or regions of Africa and nine countries of Asia [5][6][7]. In January 2019, S. frugiperda was first discovered in a cornfield of Jiangcheng City, Yunnan Province, China. Subsequently, it has invaded 29 provinces and posed a severe threat to the food security and agricultural production of China [2,4,6,8]. Quantitative real-time PCR (qRT-PCR) is an important molecular technology with high sensitivity, reliability, and specificity. It is always used in gene quantitative expression analysis and transcriptome verification [9][10][11][12]. However, many factors can lead to systematic errors and affect the accuracy of qRT-PCR results, such as the integrity and purity of total RNA, the efficiency of reverse transcription and PCR program, and pipetting [13][14][15]. Therefore, it is essential to use reference genes as the benchmark to eliminate various errors and normalize the qRT-PCR analysis [16]. Several traditional reference genes are always used as the benchmark to evaluate the expression profile of function genes, such as 18S ribosomal RNA (18S), alpha-tubulin (α-TUB), beta-actin (β-ACT), and glyceraldehyde-3phosphate dehydrogenase (GAPDH) [9,10,17,18]. These reference genes are involved in the metabolism process of insect cells. The ideal reference genes should be highly stable in various conditions [19]. So far, the reference genes have been selected in many Lepidoptera insects, including Spodoptera exigua (Hübner), Heliconius numata (Cramer), Mythimna separata (Walker), Danaus plexippus (Linnaeus), Bombyx mori (Linnaeus), Chilo suppressalis (Walker) and Helicoverpa armigera (Hübner) [20][21][22][23][24][25][26]. Furthermore, Boaventura et al. (2020) evaluate the expression stability of 10 reference genes across strains, larval instars, and gut tissues of S. frugiperda [27]. These researches proved that the reference genes show significant various expression levels under different condition treatments, and none is suitable for all experimental treatments. Accordingly, the stability of specific reference genes for various treatments should be verified before conducting qRT-PCR reactions. However, there have been only a few studies on selecting reference genes in S. frugiperda under different experimental conditions.
Here, we selected 10 candidate reference genes, 18S, AK, RPL10, RPS24, 28S, SOD, ATP, GAPDH, ACT, and a-TUB, and evaluated their stability under different treatments (including developmental stages, tissues, mating, hormones, diets, temperature) by Delta Ct method, geNorm, NormFinder, BestKeeper, and RefFinder. These findings will improve the accuracy of gene expression analysis, which can help to form a solid foundation for future research of the functional genes of S. frugiperda.

Insect Rearing
Spodoptera frugiperda (J. E. Smith) samples were obtained from Henan Academy of Agricultural Sciences, Zhengzhou, China. The larvae were fed with corn seedlings, and adults were reared with 10% sucrose solution. The rearing conditions were as follows, 26 ± 1 • C, 65% ± 5% relative humidity, and 16 h light and 8 h dark photoperiods.

Different Developmental Stages
In this study, eggs (20 mg), first-third-instar larvae (50, 20, 5 individuals, respectively), fourth-sixth instar larvae (1 individual), third-day pupae (1 male and female individual), third-day adults (1 male and female individual) were collected as samples of different developmental stages. There has a longitudinal crack with a tumor on each side of the border in the ventral side of the 9th abdominal segment of male pupae. For female pupae, there has a longitudinal crack in the ventral side with flat sides and no humps in the 8th abdominal segment.

Mating Status
To investigate the effect of mating status, mated adults of male and female (1 individual) and unmated adults of male and female (1 individual) were collected, respectively.

Diet Treatment
Spodoptera frugiperda were reared with corn seedlings, rice seedlings, wheat seedlings, and artificial feed, respectively. Then third-instar larvae (5 individuals), third-day male and female pupa (1 individual), third-day male and female adult (1 individual) were collected in each condition.

Temperature Treatment
For temperature treatment, S. frugiperda were reared at 20, 26, and 32 • C, respectively. Then third-instar larvae (5 individuals), third-day male and female pupae (1 individual), third-day male and female adult (1 individual) were collected in each temperature treatment.
The sample of each treatment was collected in 1.5 mL centrifuge tubes and rapidly frozen in liquid nitrogen, then stored at −80 • C. All treatments were set to three biological replicates.

RNA Isolation and cDNA Synthesis
According to the instruction, the total RNA of all samples was extracted by the RNAprep Pure Tissue Kit (Tiangen, Beijing, China). The RNA concentration and purity were detected with a NanoDrop spectrophotometer (MD2000C; Biofuture, UK). Subsequently, each sample of 1 µg total RNA was used for cDNA synthesis with PrimeScript RT Reagent Kit (Takara, Dalian, China).

Reference Genes Selection and Primer Design
Ten commonly used reference genes (18S, AK, RPL10, RPS24, 28S, SOD, ATP, GAPDH, ACT, and a-TUB) were selected. The primer of each gene was designed by DNAMAN V6. The purified PCR product was used as a starting template to draw the standard curve to determine the amplification efficiency of primer, and each gradient was diluted 10-fold, with 5 gradients. The melting curve had only one single peak determined that the primer was specific.

qRT-PCR Analysis
The qRT-PCR amplification was performed by C1000 Touch Thermal Cycler (Bio-Rad, Hercules, CA, USA). A total of 20 µL reaction volume containing 10 µL Fast Super EvaGreen ® Master Mix (US Everbright Inc., Suzhou, China), 1 µL cDNA of a sample, 0.5 µL of each primer (10 µM), and 8 µL ddH 2 O. The PCR reaction conditions were as follows: 95 • C for 2 min, and 40 cycles of 95 • C for 5 s, 56 • C for 30 s, 72 • C for 30 s. The melting curves of amplicons were determined by taking continuous fluorescence readings with increasing temperatures from 65 to 95 • C. Three biological repeats were set for each reaction, and three technical repeats were set for each biological repeat.

Stability Analysis
The Delta Ct method [11], geNorm [28], NormFinder [15], and BestKeeper [29] were used to evaluate the stability of each reference gene in different treatments. Based on the relative pairwise comparisons of genes within each sample, the Delta Ct method can identify the most stable reference gene. geNorm is used to calculate the M value to evaluate Insects 2021, 12, 902 4 of 13 the expression stability of reference genes. Furthermore, the geNorm can define the optimal number of reference genes by calculated pairwise variance values (V). When the value (Vn/Vn+1) is less than 0.15, it means that the optimal number is n. NormFinder is used to calculate the changes of intergroup and intragroup to determine the optimal reference gene. The BestKeeper focuses on the standard deviation and variation coefficient of Ct value, which can rank the stability of each gene. Finally, the online tool RefFinder can determine the comprehensive ranking of reference genes by integrating four algorithms.

Stability Verification of Candidate Reference Genes
Odorant-binding proteins 1 (SfruOBP1) was selected to verify the stability of reference genes. The relative expression of SfruOBP1 in tissues of S. frugiperda was calculated by the 2 -∆∆Ct method, using the most stable and unstable reference gene, respectively [30]. The significant differences were using IBM SPSS Statistics version 22 using one-way ANOVA analysis followed by Tukey's test.

qRT-PCR Analysis
Before evaluating the suitability of the reference genes, the specificity and efficiency of PCR amplification should be first confirmed. In this study, all the PCR products amplified were detected with 0.8% agarose gel, and the single band was observed in each PCR product with the expected band size. In addition, the PCR amplification for each primer pair showed a single peak in melting curves ( Figure S1). The amplification efficiency of each primer pair ranged from 90.6% to 107.1%, and all regression coefficients were higher than 0.990 (Table 1). Taken together, the primers can be used for quantitative determination. We determined the Ct value of 10 candidate reference genes for various treatments, i.e., different developmental stages ( Figure 1A We determined the Ct value of 10 candidate reference genes for various treatments, i.e., different developmental stages ( Figure 1A), various tissues ( Figure 1B), mating status ( Figure 1C), hormone treatment ( Figure 1D), diet treatment ( Figure 1E), and temperature treatment ( Figure 1F). The Ct values ranged from 11.85 (18S) to 39.12 (28S), and the Ct values of 18S were the lowest (11.85-17.25), and the 28S were the highest (25.39-39.12). The variation of Ct values in 18S was minimum (5.40), and the 28S had the maximum Ct value variation (13.73) ( Figure 1G).

Expression Stability Analysis of Candidate Reference Genes
For different developmental stages, the results of the Delta Ct method, geNorm, and NormFinder analysis indicated that SOD, RPL10, and RPS24 were the most stable gene. The BestKeeper analysis identified 18S as the most suitable reference genes ( Table 2). Following the RefFinder analysis, the overall ranking of expression stability was as follow Figure 2A). Moreover, geNorm analysis indicated that V7/V8 was less than 0.15, but based on the geNorm instructions, the top three genes were regarded as the optimal reference genes combination (SOD, RPL10, RPS24) ( Figure 3).
RPL10 was the most stable gene based on the Delta Ct method and NormFinder in various tissues. The estimation of geNorm showed that α-TUB and GAPDH were the most stable genes. However, the expression stability of ATP was highest in BestKeeper. Meanwhile, 28S was regarded as the most unstable gene based on all algorithms (Table 2). Combining four algorithms, the comprehensive ranking by RefFinder was as follow α-TUB Figure 2B). In addition, all pairwise variance values were higher than 0.15 by geNorm analysis. Thus, α-TUB, RPL10, ATP were the optimal combination of reference genes (Figure 3).
Based on the results of the Delta Ct method, geNorm, and NormFinder, AK was identified as the most stable gene in mating status; the BestKeeper analysis showed that RPL10 had the highest expression stability ( Table 2). The RefFinder calculated comprehensive ranking was as followed Figure 2C). Meanwhile, all pairwise variance values were greater than 0.15 by geNorm analysis. Therefore, the best group of reference genes for mating status were AK, RPL10, 18S (Figure 3). For hormone treatment, the Delta Ct method, NormFinder, and BestKeeper identified ATP and AK as the most suitable reference genes, and geNorm identified 18S and AK as the most suitable reference genes (Table 3). Following the RefFinder analysis, the comprehensive ranking of expression stability was as follow ATP > AK > 18S > RPL10 > RPS24 > ACT > SOD > 28S > GAPDH > α-TUB ( Figure 2D). In addition, the geNorm analysis calculated that V2/V3 was less than 0.15, indicating that two reference genes    RPL10 was the most stable gene based on the Delta Ct method and NormFinder in various tissues. The estimation of geNorm showed that α-TUB and GAPDH were the most stable genes. However, the expression stability of ATP was highest in BestKeeper. Meanwhile, 28S was regarded as the most unstable gene based on all algorithms ( Table 2). Combining four algorithms, the comprehensive ranking by RefFinder was as follow α- Figure 2B). In addition, all pairwise variance values were higher than 0.15 by geNorm analysis. Thus, α-TUB, RPL10, ATP were the optimal combination of reference genes (Figure 3).
Based on the results of the Delta Ct method, geNorm, and NormFinder, AK was identified as the most stable gene in mating status; the BestKeeper analysis showed that RPL10 had the highest expression stability ( Table 2). The RefFinder calculated comprehensive ranking was as followed AK > RPL10 > 18S > SOD > ATP > ACT > α-TUB > RPS24 > GAPDH > 28S ( Figure 2C). Meanwhile, all pairwise variance values were greater than 0.15 by geNorm analysis. Therefore, the best group of reference genes for mating status were AK, RPL10, 18S (Figure 3).
For hormone treatment, the Delta Ct method, NormFinder, and BestKeeper identified ATP and AK as the most suitable reference genes, and geNorm identified 18S and AK as the most suitable reference genes (Table 3). Following the RefFinder analysis, the comprehensive ranking of expression stability was as follow ATP > AK > 18S > RPL10 > RPS24 > ACT > SOD > 28S > GAPDH > α-TUB ( Figure 2D). In addition, the geNorm analysis calculated that V2/V3 was less than 0.15, indicating that two reference genes should be combined for normalization. So, ATP and AK were the optimal combinations of reference genes in hormone treatment (Figure 3).   For diet treatment, the Delta Ct method and geNorm indicated that the most stable gene was SOD. The expression stability of 18S was the highest based on NormFinder and BestKeeper. Meanwhile, the expression stability of 28S was lowest by four methods ( Table 3). The RefFinder calculated comprehensive ranking of the expression stability of candidate reference genes followed as 18S > RPL10 > SOD > ACT > AK > RPS24 > GAPDH > α-TUB > ATP > 28S ( Figure 2E). Moreover, the V6/V7 was smaller than 0.15 by geNorm analysis. So, the optimal combination of reference genes in diet treatment were 18S, RPL10, SOD (Figure 3).
Based on the results of the Delta Ct method, NormFinder, and BestKeeper, RPL10 was identified as the most stable gene in temperature treatment, whereas geNorm identified RPS24 and SOD as the most suitable reference genes ( Table 3). The stability comprehensive order of these genes in RefFinder analysis was listed as following RPL10 > 18S > RPS24 > ACT > SOD > AK > α-TUB > GAPDH > ATP > 28S ( Figure 2F). Furthermore, geNorm analysis indicated that V7/V8 is less than 0.15, RPL10, 18S, and RPS24 were regarded as the optimal combination of reference genes for temperature treatment (Figure 3).
For all samples, RPL10 was regarded as the most stable gene based on the Delta Ct method, geNorm, and NormFinder. The expression stability of 18S was the highest in Best-Keeper analysis. Meanwhile, the stability of 28S was the worst in all algorithms (Table 3). Combining four algorithms, the comprehensive order by RefFinder was as follow RPL10 > 18S > SOD > ACT > AK > GAPDH > α-TUB > ATP > RPS24 > 28S ( Figure 2G). In addition, the V7/V8 was smaller than 0.15 by geNorm analysis. So, the optimal combination of reference genes in all samples were RPL10, 18S, SOD (Figure 3).

Verification of Candidate Reference Genes
We selected SfruOBP1 as the target gene, and four candidate genes (α-TUB, RPL10, ATP, 28S) were used as reference genes to determine the expression level of the target gene in different tissues of S. frugiperda. When the more stable genes (α-TUB, RPL10, ATP) were used as the reference genes, SfruOBP1 was significantly expressed in the head of females ( Figure 4A-C). However, when the unstable gene (28S) was used as the reference gene, the expression level of SfruOBP1 in the head of males was significantly higher than females ( Figure 4D). ATP, 28S) were used as reference genes to determine the expression level of the target gene in different tissues of S. frugiperda. When the more stable genes (α-TUB, RPL10, ATP) were used as the reference genes, SfruOBP1 was significantly expressed in the head of females ( Figure 4A-C). However, when the unstable gene (28S) was used as the reference gene, the expression level of SfruOBP1 in the head of males was significantly higher than females ( Figure 4D).

Discussion
The qRT-PCR is a reliable technique in gene expression analysis, with high sensitivity and specificity, and the qRT-PCR data must be normalized by suitable reference genes to avoid expression differences among samples [31]. Recently, there have many reports on reference genes selection of various insects under different abiotic and biotic conditions [32][33][34][35][36], including S. frugiperda [27]. In these studies, the expression levels of conventional reference genes in different insects are quite different, and none of them with similar expression levels under all conditions, which indicates that there is no absolute universality among homologous reference genes. Thus, screening reference genes are essential for quantitative research under certain conditions.
In this study, the Delta Ct method, geNorm, NormFinder, and BestKeeper were used to evaluate the expression stability of candidate reference genes. The results indicated that the ranks of candidate reference genes were significantly different by various algorithms. For example, the expression stability of RPL10 was the highest in various tissues of S. frugiperda based on the Delta Ct method and NormFinder. Moreover, the geNorm and BestKeeper indicated that the most stable gene was GAPDH and ATP, respectively. In addition, RPL10 was regarded as the most stable gene based on the Delta Ct method, NormFinder, and BestKeeper under hormone treatment, and the expression stability of 18S was the highest in geNorm analysis. Similarly, following the Delta Ct method, geNorm, and NormFinder analysis, Xu et al. (2017) find that the expression stability of EF1 is the highest of Chilo suppressalis (Walker) under temperature stress, and RPS11 is the lowest, but RPS11 is the most stable gene, and EF1 is an unstable gene in BestKeeper analysis. These results also occur in other insects, as each program is based on its own unique algorithm [24,37]. Therefore, in order to eliminate the error of the algorithm of the four methods, the online program RefFinder was used to provide a comprehensive ranking. The results showed that SOD for different developmental stages, α-TUB for various tissues, AK for mating treatment, 18S for hormone treatment, ATP for diets treatment, RPL10 for temperature treatment were the most stable reference genes of S. frugiperda. These results also further verified that the expression stability of the reference gene varied under different treatments.
Ribosomal protein (RP) is the component of ribosomes. It plays a vital role in protein biosynthesis in all biological cells. It also participates in cell growth regulation, cell differentiation, and DNA repair [38,39]. In this paper, RPL10 was the most stable gene in temperature treatment and all test samples of S. frugiperda. Meanwhile, RPL10 also showed a relatively high stable ranking in other treatments of S. frugiperda. Similarly, Boaventura et al. evaluated several candidate reference genes for normalization of gene expression data across strains, larval instars, and gut tissues of S. frugiperda and identified that RPL10 was a stable reference gene [27]. Previous research also indicates that ribosomal protein genes, as commonly used reference genes, have a wide range of applications in many insects. For example, RPL27 and RPL32 in different developmental stages and tissues of Apolygus lucorum (Meyer-Dur) [33], RPL12 in starvation treatment of Anthonomus eugenii (Cano) [9], and RPS20 in temperature treatment of Diaphorina citri (Kuwayama) [40]. However, the stability of ribosomal protein genes is not universal in all insects [24,34].
With the in-depth study of the stability of reference genes, many researchers have proposed that using multiple reference genes can improve the accuracy of qRT-PCR and remove biased normalization [15,41,42]. The geNorm analysis not only evaluates the stability of reference genes but also calculates the best combination of the reference gene under certain conditions. Based on the geNorm instructions, when pairwise variance value is less than 0.15, it means that the number of optimal combinations is n; when the pairwise variation values are greater than 0.15, we can select two or three most stable reference genes as combination according to the trend of pairwise variation value. In this paper, the pairwise variance value V2/V3 of hormone treatment was less than 0.15, so we needed two reference genes to analyze the gene expression. Although the V3/V4 of other treatments were more than 0.15, we propose that three of the best reference genes should be used to obtain accurate and reliable results under these conditions. Many studies have shown that more reference genes are needed for accurate normalization when the sample size becomes larger. However, the accuracy of target gene expression is decreased when the fourth reference gene is introduced [41,42]. It followed that the threshold of pairwise variance values (V < 0.15) is not absolute.
Odorant-binding protein is one of the most important proteins in the olfactory system of insects [43]. In this study, we determined the expression levels of SfruOBP1 in various tissues to verify the stability of reference genes. The results indicated that the expression patterns of SfruOBP1 were similar when the most stable and recommended combination of reference genes in various tissues were used as the benchmark, respectively; it was significantly expressed in the head of the female. However, the expression pattern of SfruOBP1 was different from that mentioned above when the least stable reference gene was selected as the benchmark, and SfruOBP1 was significantly expressed in the head of the male. Therefore, screening suitable reference genes are essential to improve the accuracy of function gene expression.

Conclusions
In conclusion, the stability of 10 candidate reference genes was analyzed by five reliable algorithms under different experimental conditions. The results showed that the optimal combination of most stable reference genes was SOD, RPL10, and RPS24 for developmental stages; α-TUB, RPL10, and ATP for various tissues; AK, RPL10, and 18S for mating status; 18S and AK under hormone treatment; 18S, RPL10, and SOD under diet treatment; RPL10, 18S, and RPS24 under temperature treatment. This study confirmed recent data on some reference genes and provided an evaluation of a number of additional reference genes of S. frugiperda under various conditions.