Evaluation of reference genes for quantitative expression analysis in Mylabris sibirica (Coleoptera, Meloidae)

Mylabris sibirica is a hypermetamorphic insect whose adults feed on oilseed rape. However, due to a shortage of effective and appropriate endogenous references, studies on molecular functional genes in Mylabris sibirica, have been tremendously limited. In this study, ten internal reference genes (ACT, ARF1, AK, EF1α, GAPDH, α-TUB, RPL6, RPL13, RPS3 and RPS18) were tested and assessed under four selected treatments including adult ages, adult tissues, temperatures, and sex by RT-qPCR based on five methods (Ct value, geNorm, NormFinder, BestKeeper and RefFinder). Our findings showed that RPL6 and RPL13 were the most optimal internal reference gene combination for gene expression during various adult ages and under diverse temperatures; The combination of RPL6 and RPS18 was recommended to test gene transcription levels under different adult tissues. AK and RPL6 were the best reference genes in male and female adults. RPL6 and RPL13 were the most appropriate reference gene pair to estimate gene expression levels under four different tested backgrounds. The relative transcript levels of a uridine diphosphate (UDP)-N-acetylglucosamine-pyrophosphorylase (MsUAP), varied greatly according to normalization with the two most- and least-suited reference genes. This study will lay the basis for further molecular physiology and biochemistry studies in M. sibirica, such as development, reproduction, sex differentiation, cold and heat resistance.


Introduction
Quantitative real-time fluorescent polymerase chain reaction (qRT-PCR) is a crucial method to measure target gene transcripts (Derveaux et al., 2010;Jozefczuk and Adjaye, 2011) and microbial abundance (Ali et al., 2018a;Ali et al., 2018b;Ali et al., 2018c;Ali et al., 2019) due to its high specificity, sensitivity and accuracy.Its precision is affected by numerous biological and technical factors, such as RNA purity, PCR efficiency, inappropriate data and statistical analyses (Valasek and Repa, 2005;Bustin et al., 2009).
When qRT-PCR is applied to test transcripts, it is indispensable to normalize to improve the quantitative results by combining relatively stable reference genes (Bustin et al., 2009).If poor internal genes are used, the quantitates of nucleic acid will be inaccurate (Yang et al., 2014).Hence, Suitable internal references should be verified under diverse backgrounds, including developing stages, tissues and hosts (Andersen et al., 2004;Zhang et al., 2022;Shen et al., 2022).
As usual, reference genes are housekeeping genes (HKGs), whose expression levels are stable during different physiological states or in different cells (Zhang et al., 2022).At present, the 10 most generally used references contain actin (ACT), glyceraldehyde-3-phosphate dehydrogenase (GAPDH), ribosomal protein (RPL and RPS), TATA binding protein (TATA), heat shock protein (HSP), elongation factor 1α (EF1α), tubulin (TUB), 18S ribosomal RNA (18S) and succinate dehydrogenase complex subunit A (SDHA) (Lü et al., 2018).In addition, several novel methods including bioinformatic transcripome analysis have been used for reference gene screening (Yu et al., 2020).In previously researches, one to five analysis tools were applied to evaluate gene expression stability (Lü et al., 2018).Currently, the RefFinder is the only web-based tool available for evaluating candidate housekeeping genes by integrating the four available computational programs (geNorm, NormFinder, BestKeeper, and ΔCt method) into a web-based tool (Xie et al., 2023).
Nevertheless, there is no one stable internal gene under numerous tested treatments, including developmental stages, tissues and hosts (Andersen et al., 2004;Zhang et al., 2022;Shen et al., 2022;Shen et al., 2023).To evaluate accurate gene expression levels, each candidate reference gene under diverse tested conditions must be validated (Xu et al., 2021).
Mylabris sibirica (Coleoptera, Meloidae) is a hypermetamorphic pest that leads to significant losses in oilseed rape (Brassica napus) production (Tian Yubo et al., 2021).Its adults mainly feed on flowers of B. napus (Tian Yubo et al., 2021).Until now, studies on M. sibirica have concentrated on classification (ČernÝ and Vrabec, 2019;Pan and Ren, 2020), phylogenetics (Bologna et al., 2005) and medical value (Tian Yubo et al., 2021), however, little is known on gene functions.Gene expression analyses are essential to study molecular mechanisms of development, reproduction and physiology.Nevertheless, no study evaluating M. sibirica gene expression analysis has been reported.In order to further manage M. sibirica based on novel target genes, screening for optimal reference genes is imperative.

Experimental conditions
The recommended reference genes

FIGURE 3
Stabilities of the ten house-keeping genes in Mylabris sibirica among various adult tissues.Foregut, midgut, hindgut, epidermis, trachea and antenna were dissected from the newly-emerged adults.The expression stability rankings are determined by geNorm, NormFinder and BestKeeper.
the stability of gene expression during various adult ages, among various adult tissues, under diverse temperatures, and in females and males in M. sibirica.Finally, uridine diphosphate (UDP)-N-acetylglucosamine-pyrophosphorylase (UAP) was used as the target gene to verify our findings.Our results will offer the reference foundation for further molecular mechanisms in M. sibirica.

Methods and materials
Insect collection   Stabilities of the ten house-keeping genes in Mylabris sibirica in males and females.The male and female adults were collected at day 5 after emerging.The expression stability rankings are determined by geNorm, NormFinder and BestKeeper.

Specimens through diverse adult ages
The newly emerged adults were grouped as males and females and kept in different net cages (30 cm × 30 cm × 40 cm, handmade).For each biological replicate, a total of five adults (three males and two females) were collected every 24 h for a continuous period of 72 h.Three biological replicates were conducted.

Samples in different adult tissues
Epidermis, foregut, midgut, hindgut, trachea and antenna were dissected from the newly-emerged adults.Ten adults (5 males and 5 females) were dissected for one biological replicate.Three biological replicates were prepared.

Specimens under various temperatures
The newly-emerged adults were reared under five temperatures (4 °C, 15 °C, 20 °C, 26 °C and 35 °C) for 8 h.Five adults (3 males and 2 females) were collected for one biological replicate.Three biological replicates were prepared.

Collections in males and females
The newly-emerged male and female adults were kept in an insectary at 26 °C for 5 days.Five males and females were collected respectively for one biological replicate.Three biological replicates were prepared.

Total RNA extraction and cDNA synthesis
Total RNAs were extracted by TRIzol reagent (YiFeiXue Tech, Nanjing, China) following the manufacturer's instructions.The RNA integrity was assessed by electrophoresis on a 1.5% agarose gel.The purity and concentration of the total RNA samples were determined using a NanoDrop ND-1000 spectrophotometer (Nanodrop Technologies, Rockland, DE, United States), with 260/280 and 260/ 230 ratios ranging from 1.9 to 2.1.The HiScript III RT SuperMix (containing gDNA wiper, Vazyme Biotech Co.,Ltd., Nanjing, China) was used to synthesize cDNA.The reaction mixtures were incubated at 37 °C for 15 min, followed by 85 °C for 5 s.The resulting cDNA samples were diluted 5-fold for the pursuant PCR and RT-qPCR.
The primers of Reverse transcriptase PCR (RT-PCR) were designed by Primer Premier 5.0 and performed according to the previously described method (Shen et al., 2022).The resultant sequences were submitted to GenBank.The accession numbers were located in Supplementary Table S1.

Quantitative real-time PCR (qRT-PCR)
The primers of qRT-PCR were designed by NCBI Primer-BLAST, and were listed in Table 1.The qRT-PCR reactions were conducted using ChamQ Universal SYBR qPCR Master Mix (Vazyme Biotech Co., Ltd.) and CFX96 Real-Time System (Bio-Rad Laboratories, Hercules, CA, United States).The reaction mixture consisted of a 20 μL reaction volume including 1 μL cDNA template, 0.4 μL of forward primer (10 μM), 0.4 μL of reverse primer (10 μM), 10 μL of 2× ChamQ Universal SYBR qPCR Master Mix and 8.2 μL of RNase Free water.The qRT-PCR reaction conditions consisted of the following steps: an initial step of 95 °C for 30 s, followed by 40 cycles of denaturation at 95 °C for 5 s, and subsequent annealing at 60 °C for 34 s, followed by one cycle of 95 °C for 15 s, 60 °C for 60 s and 95 °C for 1 s.The specificity of PCR amplicons were assessed by a melting curve analysis, analyzing by the BioRadCFXManager and gel electrophoresis.All experiments were repeated in triplicate.Amplification efficiencies (E) were calculated using a 5-fold dilution series of template via the following equation: E = (10 [−1/slope] − 1) × 100%.

Stability determination of candidate reference genes
Uridine diphosphate (UDP)-N-acetylglucosaminepyrophosphorylase (UAP) of M. sibirica was used to verify the stability of candidate reference genes (GenBank: OR838722).The primer sequence of the target gene was as follows: Forward: ATTATTGATGGCCGGTGGTC Reverse: ACCATTTAAACCGGTCTTTTGTT Based on the stability (RPL6 and RPS18) and instability (AK and GAPDH) of primary reference genes, the relative levels of MsUAP in different adult tissues were computed by the 2 −ΔΔCT method and from three replicates.One-way analysis of variance was used to assess significance in UAP expression levels among various adult tissues (SPSS, Chicago, IL, United States).

Data processing
The raw Ct values were obtained using the BioRadCFXManager.The stability of candidate HKGs were measured using the ΔCt method (Silver et al., 2006), geNorm (Vandesompele et al., 2002), Normfinder (Andersen et al., 2004) and BestKeeper (Pfaffl et al., 2004).Furthermore, the optimal number of reference genes for gene expression normalization was determined by pairwise variation (V n/ n+1 ) using the GeNorm program Typically, a V n/n+1 value below the threshold of 0.15 indicates that the starting n reference genes are enough for normalizing target gene expression.Lastly, the overall ranking of each experimental background was assessed based on RefFinder (Xie et al., 2012;Xie et al., 2023).
The products obtained from qRT-PCR were validated through sequencing.The primer specificity for qRT-PCR was confirmed through melting curve analysis.As expected, slopes of all primer pairs were less than −3.0, and regression coefficient (R2) and efficacy values ranged from 0.992-1 and 90.94%-99.25%,respectively (Table 1).Our findings demonstrated that efficiency met the required standards for traditional qRT-PCR (Bustin et al., 2009).

Stability of the ten HKGs
Samples were collected from three adult ages (the first to third days of the newly-emerged adults), six adult tissues (foregut, midgut, hindgut, epidermis, trachea and tentacle), five temperature treatments (4 °C, 15 °C, 20 °C, 26 °C and 35 °C) and both sexes (5 days after emerging).We discovered that the mRNA levels of ten HKGs were abundant at various adult ages, during diverse adult tissues, under different temperatures and both sexes using qRT-PCR,.
The all cycle threshold values (Ct) under diverse tested backgrounds were presented (Figure 1).The ΔCt method Stabilities of the ten house-keeping genes in Mylabris sibirica in different samples.The stability of the reference genes as calculated by the Geomean method of RefFinder.A lower Geomean of ranking value denotes more stable expression.(A) Adult ages, (B) adult tissues, (C) temperature, (D) sex, (E) all samples.
assesses the genes stability based on genes average of STDEV (Silver et al., 2006).During diverse adult ages, the expression fluctuations of RPL6 and AK were smaller, whereas ACT and GAPDH were higher (Figure 1A).Under diverse tissues, except for AK, the expression variations were small in ten HKGs (Figure 1B).Under various temperatures, the expression difference of RPL13 and ARF1 was smaller (Figure 1C).In males and females, the expression variations of AK and EF1α were small in ten HKGs (Figure 1D).

Stability of the ten HKGs under adult ages
The geNorm statistical analysis assesses the gene stability by M-values (the average expression stability) and V-values (pairwise variation).These results showed that EF1α and RPL13 were the most stable internal genes during various adult ages, with M-values below 0.25.In contrast, the most unstable genes were ACT and GAPDH (Figure 2A; Table 2).Pairwise variation analysis showed that all values were below 0.15, demonstrating that two reference genes were needed for the gene expression analysis under adult ages (Figure 2B).
Based on the NormFinder algorithm, the stable rankings of ten internal genes from the most stable to the least were RPL6, AK, α-TUB, RPS3, RPL13, RPS18, EF1α, ARF1, RPL27, GAPDH and ACT (Figure 2C; Table 2).Except for GAPDH and ACT, the p values of other genes were had less than 1 (Figure 2C; Table 2).

Stability of the ten HKGs during different adult tissues
According to the geNorm data, the reference genes rankings ranging from the most stable to the least stable were RPS3, RPS18, RPL13, RPL6, EF1α, ARF1, ACT, α-TUB, GAPDH and AK (Figure 3A; Table 2).Pairwise variation analysis showed that all values were less than 0.15, indicating that two reference genes from different tissues were needed for the gene expression analysis (Figure 3B).
The BestKeeper analysis showed that RPL6, were the most stable, with the SD values of 0.604 (Figure 3D; Table 2).Again, except for AK, the SD values of other reference genes were below 1.0 (Figure 3D; Table 2).

Expression stability of the ten HKGs under various temperature treaments
The geNorm analysis showed that RPL6, RPL13, ARF1 and EF1α were the most stable internal genes under different temperatures, whose M-values were 0.233, 0.233, 0.259 and 0.283, respectively (Figure 4A; Table 2).Moreover, pairwise variation analysis showed that all values were less than 0.15, indicating that two different internal genes are needed for testing gene expression during various temperature treaments (Figure 4B).The NormFinder data uncovered that the steady rankings were RPL13, ARF1, RPS3, RPL6, α-TUB, GAPDH, EF1α, RPS18, AK and ACT (Figure 4C; Table 2).

Stability of the ten HKGs in males and females
Based on the geNorm algorithm results, RPL6, RPL13, RPS18, EF1α and ARF1 were the most stable, with the M-values less than 0.5 (Figure 5A; Table 2).In addition, pairwise variation data manifested that the V2/3 to V5/6 values were less than 0.15, indicating that two different internal genes are suitable for evaluating gene expression levels in males and females (Figure 5B).
The BestKeeper result demonstrated that AK, EF1α, RPL6 and RPL13 were the most optimal genes, whose SD values were 0.109, 0.345, 0.469 and 0.484, respectively (Figure 4D; Table 2).Again, values of these genes were less than 0.5, showing their similar stabilities.

Validation of the selected reference genes
To assess the stability of the candidate reference genes, the relative expression level of UAP was assessed in the epidermis, foregut, midgut, hindgut, trachea and antenna.The following reference genes were used to normalize: RPL6 + RPS18 (the most stable reference gene), and AK and GAPDH (the least stable reference gene).The highest accumulated mRNA level of MsUAP was detected in the midgut, hindgut and trachea, followed by those in the epidermis, the lowest level was found in the foregut and antenna.However, AK and GAPDH was used as reference genes, and MsUAP was highly expressed in the trachea, lowly expressed in the epidermis (Figure 7).

Discussion
Quantitative real-time PCR (qRT-PCR) is a popular method for assessing gene expression (Derveaux et al., 2010;Jozefczuk and Adjaye, 2011) and microbial abundance determination (Ali et al., 2018a;Ali et al., 2018b;Ali et al., 2018c;Ali et al., 2019).Optimal reference genes are critical greatly in eliminating heterogeneity in diverse datasets and improving the quantitative results (Bustin et al., 2009).An ideal internal gene are abundantly expressed under diverse experimental treatments (Derveaux et al., 2010).However, not all internal genes remains suitable in various species (Zhang et al., 2022).Therefore, selecting suitable references must be conducted before qRT-PCR.Studies on internal genes have been reported in many insect species (Shakeel et al., 2018), such as Aphidoletes aphidimyza (Shen et al., 2023), Nilaparvata lugens (Zhang et al., 2022), Phthorimaea operculella (Shen et al., 2022), H. vigintioctomaculata (Zhang et al., 2022) and P. striolata (Guo et al., 2023).As an oilseed rape pest, evaluating optimal reference genes are helpful in studying molecular mechanisms in M. sibirica.As we all knows, the current research is the first report on reference genes assessment in M. sibirica.
In male and female adults, the most reliable reference genes were AK and RPL6 (Figures 1, 5, 6; Table 2).As we discussed in the previous paragraph, RPL6 was expressed abundantly in males and females.Arginine kinase (AK) has a prominent function in invertebrate energy metabolism, which catalyzes the reversible phosphorylation of l-arginine to make up phosphoarginine (Shi et al., 2012).Similar to our results, the arginine kinase gene has been selected as the most suitable gene for expression assessment in Xylosandrus germanus (Patwa et al., 2021), Spodoptera frugiperda (Han et al., 2021) and Spodoptera litura (Lu et al., 2013).
Moreover, The BestKeeper results demonstrated the SD values of ACT, α-TUB and GAPDH were above 1.0 (Table 2), suggesting that the three reference genes were unsuitable as internal genes for RT-qPCR.
Actin (ACT) is extremely abundant in eukaryotes, and has an vital function in cellular activities, including cell motility and the regulation of transcription cell secretion (Dominguez and Holmes, 2011).At present, many publictions have verified that the expression level of ACT is less steady in various insects, such as P. operculella (RPL13) (Shen et al., 2022), H. vigintioctopunctata (RPL13 and RPS18) (Lü et al., 2018) and Hippodamia convergens (Yang et al., 2016).
Microtubule α-tubulin (α-TUB), interacts with many microtubule-associated proteins to conduct a variety of cellular functions, such as intracellular transport and cell division (Hammond et al., 2008).Currently, many studies indicated that the stability of GAPDH expression was greatly unsteady in P. operculella (Shen et al., 2022), H. vigintioctomaculata (Zhang et al., 2022), N. lugens (RPL5, RPS8 and RPL14) (Zhang et al., 2022) and A. aphidimyza (Shen et al., 2023), In brief, our results recommend RPL6 and RPL13 as the most stable internal gene pair under four tested backgrounds (Figure 6; Table 3).To further validate the reference genes in M. sibirica, we assessed the relative expression level of UAP in different adult tissues.Our findings revealed that UAP expression pattern was inconsistent in the different adult tissues when normalized to the two best-and least-stable reference genes (Figure 7).These results showed that the unreasonable use of reference genes may lead to inaccurate results for target genes.Therefore, it is crucial to choose and validate the best reference genes to ensure the accuracy of gene expression.Our study would provide a foundation for further gene molecular functions and microbial abundance determination in M. sibirica.

FIGURE 1
FIGURE 1 Expression levels of ten house-keeping genes in Mylabris sibirica.The mean C t values for 10 candidate reference genes are shown in four different experiments: (A) adult ages, (B) adult tissues, (C) temperature, (D) sex.Mean Ct values for the eight candidate reference genes are presented in box plots, where each box indicates the 25th and 75th percentiles, and the line across the box represents the median.Abbreviation: ACT, actin; ARF1, ADP ribosylation factor 1; AK, arginine kinase; α-TUB, α-tubulin; GAPDH, glyceraldehyde-3-phosphate dehydrogenase; EF1α, elongation factor 1α; RPL6, RPL13, RPS3 and RPS18, ribosomal protein.The abbreviations are exactly the same as Figure 2 through Figure 6.

FIGURE 2
FIGURE 2Stability of the ten house-keeping genes in Mylabris sibirica during diverse adult ages.All stages of Mylabris sibirica adults were sampled (collected on the first to third days of the newly-emerged adults).The expression stability rankings are determined by geNorm, NormFinder and BestKeeper. et

FIGURE 4
FIGURE 4 Stability of the ten house-keeping genes in Mylabris sibirica under different temperatures.The newly-emerged adults reared under five temperatures (4 °C, 15 °C, 20 °C, 26 °C and 35 °C) were collected.The expression stability rankings are determined by geNorm, NormFinder and BestKeeper.

FIGURE 7
FIGURE 7Relative gene expression of UAP in different adult tissues of Mylabris sibirica.The relative gene expression level of UAP in the foregut, midgut, hindgut, epidermis, trachea and antenna was normalized to the best stable [(A) RPL6 and RPS18] and least stable [(B) AK and GAPDH] reference genes, respectively.The values are means +SE.Different letters indicate significant differences in gene expression among different tissues (p < 0.05).

TABLE 1
Primers of 10 candidate house-keeping genes used in qRT-PCR.

TABLE 2
Expression stability of the candidate reference genes under different experimental conditions.

TABLE 2 (
Continued) Expression stability of the candidate reference genes under different experimental conditions.
a Candidate reference gene.

TABLE 3 A
list of the recommended reference genes in M. sibirica for different experimental conditions.