Selection of Reference Genes for Normalization of Real-Time PCR Data in Calliptamus italicus (Orthoptera: Acrididae) Under Different Temperature Conditions.

Global warming has dominated worldwide climate change trends, and adaptability to high temperatures is the main factor underlying the spread of the pest Calliptamus italicus in Xinjiang Province, China. However, knowledge about the molecular mechanisms responsible for this adaptability and other related biological properties of C. italicus remain relatively unclear. Real-time quantitative polymerase chain reaction (RT-qPCR) is a key tool for gene expression analysis associated with various biological processes. Reference genes are necessary for normalizing gene expression levels across samples taken from specific experimental conditions. In this study, transcript level of five genes (GAPDH, 18S, TUB, ACT, and EF1α), commonly used as reference genes, were evaluated under nine different temperatures (27, 30, 33, 36, 39, 42, 45, 48, and 51°C) to assess their expression stability and further select the most suitable to be used on normalization of target gene expression data. Gene expression profiles were analyzed using geNorm, NormFinder, and BestKeeper software packages. The combined results demonstrated that the best-ranked reference genes for C. italicus are EF1α, GAPDH, and ACT under different thermal stress conditions. This is the first study that assesses gene expression analysis across a range of temperatures to select the most appropriate reference genes for RT-qPCR data normalization in C. italicus. These results should assist target gene expression analysis associated with heat stress in C. italicus.

Calliptamus italicus is a pest that occasionally has a substantial impact on crops (Blanchet et al. 2010). Global air temperatures have risen by 0.6°C in the past 100 yr, and this trend is expected to persist into the future (Walther et al. 2002, Rogelj et al. 2016). Notably, insects are highly sensitivity to high environmental temperatures (Wang et al. 2014, Nishide et al. 2015, Warren et al. 2018), and community responses to global warming will likely retain strong physiological signals (Hamblin et al. 2017). Calliptamus italicus is found in Xinjiang Province, China, especially in Chapchal County, Yili State, and it is one of the predominant insect species throughout the province. Calliptamus italicus has caused substantial damage to grassland ecosystems (Xu et al. 2019). It is possible that there should be some physiological change in C. italicus as a consequence of global warming. Therefore, elucidating the mechanisms of physiological and stress adaptation to thermal conditions for the species is important (Wang et al. 2014), including the expression of key genes associated with adaptation to thermal conditions.

Gene expression studies yield important insights on physiological responses and assist in the interpretation of other research. Real-time quantitative polymerase chain reaction (RT-qPCR) methods have become a most widely applied technology for studying gene expression levels (Livak and Schmittgen 2001). These methods are relatively accurate, fast, easy-to-use, and sensitive, thereby enabling the monitoring of gene transcription (Heid et al. 1996, Wang et al. 2014). However, the means of executing the precise quantification of template nucleic acid is a central problem in the analysis of gene transcription. There are numerous critical components in the workflow that need to be accounted for in order to reach biologically meaningful and trustworthy conclusions (Derveaux et al. 2010, Aggarwal et al. 2018). Generally, the most common approach to normalizing data from gene expression experiments is to use reference genes as internal controls (Ye et al. 2018). Accordingly, it is crucial that an appropriate reference gene is chosen (Radonic et al. 2004, Zhao et al. 2019). The ideal reference gene should not be regulated or influenced by experimental conditions. However, there has not yet been discovered a reference gene with universally stable expression across all experimental conditions (Yang et al. 2014b). Usually, reference genes, such as glyceralde-hyde-3-phosphate dehydrogenase (GAPDH), the 18S ribosomal RNA (18S) subunit, α -tubulin (TUB), and actin (ACT) (Nicot et al. 2005, Zhao et al. 2019), are selected for normalization of gene expression.

Although, a previous study has already demonstrated the most suitable reference genes across different developmental stages for gene profile studies in Locusta migratoria, e.g., EF1α, Hsp70, and RPL32 (Yang et al. 2014a), there have been no systematic screenings of reference genes across different temperatures for C. italicus. In this study, we selected five candidate reference genes (GAPDH, 18S, TUB, ACT, and EF1α) and analyzed their expression stability across nine different temperatures. The geNorm (Vandesompele et al. 2002), NormFinder (Andersen et al. 2004), and BestKeeper (Pfaffl et al. 2004) statistical software programs were used to assess which gene among GAPDH, 18S, TUB, ACT, and EF1α were the most stable and therefore represents the best choice for normalizing data in gene expression studies of C. italicus (Orthoptera: Acrididae).

Materials and Methods

Experimental Insects and Exposure to Temperature Stress

Adult C. italicus individuals were collected in the field (Chapchal County, Yili State, Xinjiang Province, China) in July 2018 and brought to the laboratory, where they were reared without food for one day at 27°C (the standard rearing temperature for the species) under a 10-h photoperiod cycle and prepared for the temperature stress assay.

Healthy C. italicus were randomly divided into nine experimental groups, each containing 30 individuals. Nine artificial climate boxes were used to control temperature, and C. italicus were exposed to nine different temperatures (27, 30, 33, 36, 39, 42, 45, 48, and 51°C; Li et al. 2015) under light. After a 2-h exposure, live C. italicus individuals were placed into liquid nitrogen for 2 h and stored at −80°C until RNA extraction.

Three biological replicates and three technical replicates per temperature were used in this study, and there were three C. italicus individuals in each biological replicate.

Total RNA Extraction and cDNA Synthesis

Total RNA was isolated from 30 mg of three C. italicus individual tissue homogenates by using liquid nitrogen in a mortar and pestle. The RNeasy mini kit (Qiagen, Hilden, Germany) was used to extract total RNA according to the manufacturer’s protocol, and RNA was subsequently purified with the RNase-Free DNase Set (Qiagen) to remove any remaining genomic DNA. RNA concentration and purity were measured using a NanoDrop2000 spectrophotometer (ThermoFisher Scientific, Waltham, MA), and RNA quality was assessed by 1% agarose gel electrophoresis. Then, 1 μg of RNA was reverse transcribed in a 30-µl volume reaction using oligo d(T)15 primer (Takara, Kusatsu, Japan) and M-MLV reverse transcriptase (Takara, Kusatsu, Japan) following the manufacturer’s protocol. Reverse transcription products were stored at −20°C.

Selection of Candidate Reference Genes and Primer Design

All five candidate reference genes were selected from published literature on locust species (Van Hiel et al. 2009, Zhao et al. 2012, Yang et al. 2014a) and were retrieved from GenBank. The primers are described in Table 1. Primers were designed using PrimerSelect software (DNASTAR, Inc., Madison, WI) based on the GenBank sequences. Primer specificity was verified by checking the PCR products via 1% agarose gel electrophoresis (Supp File 1 [online only]). Melt curves were generated based on temperatures ranging from 65 to 95°C to check the specificity of the reactions.

Table 1.

List of primers and their characteristics

Gene	Full name	Primers sequence（5′–3′）	Amplication length (bp)	Regression coefficients (R2)	Eficiencies (E%)	Mean Cq ± SD	Coefficient of variation (CV%)	Function	GenBank accession numbers
GAPDH	Glyceraldehyde-3-phosphate dehydrogenase	TGAAATTGTTGAGGGATTGATGA	178	0.99	95	17.83 ± 0.86	5.11%	Oxydoreductase in glycolysis	MN421945
		CACTGGAACTCTGAAAGCCAT
								and gluconeogenesis
EF1	Elongation factor 1 alpha	GTGGGCCGAGTAGAAACAGG	198	0.99	116	20.60 ± 0.63	3.25%	Translation eukaryotic factor	MN421944
		TGAATCACCAGCAACATAACCAC
TUB	Alpha-tubulin	CGAGCCATACAATTCCATCCTTAC	162	1.00	105	30.56 ± 2.69	9.28%	Cytoskeletal structural protein	MN421943
		GAAACTATCTGGCCAATCAACCTG
ACT	β-Actin	CGAAACCTTTAATACCCCAG	102	0.99	108	22.28 ± 0.46	2.20%	Cytoskeletal structural protein	MN421941
		CCATCACCAGAATCCAACAC
18S	18S ribosomal RNA	ATGCAAACAGAGTCCCGACCAGA	154	0.99	105	15.68 ± 0.48	3.24%	Structural RNA for the small	MN421942
		CCTGGTGGTGCC CTTCCGTCAA
								Component of eukaryotic
								Cytoplasmic ribosomes.

RT-qPCR Assay

LightCycler 480 SYBR Green I Master (Roche, Basel, Switzerland) on the LightCycler 96 System (Roche) was used for RT-qPCR. Each 20-µl reaction consisted of 10 µl of LightCycler 480 SYBR Green I Master (2× concentration), 1 µl of cDNA template, 0.3 µl of 10 μmol/liter forward and reverse primers each, and 8.4 µl of nuclease-free water. Reactions were conducted at 95°C for 10 s as an initial denaturation step, followed by 40 cycles of 95°C for 5 s, 60°C for 30 s, and 72°C for 15 s. Three technical replicates were analyzed for each biological sample.

The primers’ regression coefficients (R2) and efficiencies of amplification (E) were determined based on the standard curves, which were constructed by amplifications of a series of six 10-fold dilutions of cDNA, and E was calculated for each primer pair by determining the slopes of the standard curves following the equation E(100%) = (10(−1/slope) − 1) × 100. The performance of each primer set was assessed according to E and R2 values (Table 1). Primer specificity was verified using a melt curve analysis. The coefficient of variation (CV) was calculated following the equation CV (%) = (SD/mean) × 100% (Radonic et al. 2004).

Data Analysis of Reference Gene Expression Stability

The expression stabilities of the five reference genes across the different temperature treatments were evaluated using three algorithms: geNorm (Vandesompele et al. 2002), NormFinder (Andersen et al. 2004), and BestKeeper (Pfaffl et al. 2004). For the geNorm and NormFinder algorithms for analyses, the raw Cq values were converted into relative quantities using the formula 2−ΔCq (ΔCq = [each sample Cq value] − [the lowest Cq value]). For the BestKeeper analysis, the average Cq values from each sample were directly used (Ye et al. 2018, Zhao et al. 2019).

NormFinder

NormFinder evaluates the expression stability of candidate reference genes on the basis of intra- and inter-group comparisons for each reference gene. The high-expression stability of this gene is reflected in a low stability value (SV) (Andersen et al. 2004).

geNorm

Similarly, geNorm calculates the expression stability of a candidate reference gene according to stability value (M). Lower M values represent stable candidate reference genes, while higher values reflect less stable genes. In addition, geNorm also evaluates pairwise variation values (V), which determines the lowest number of reference genes needed for accurate normalization (Vandesompele et al. 2002, St-Pierre et al. 2017).

BestKeeper

BestKeeper determines the best reference genes based on standard deviation (SD), Pearson correlation coefficient (r), and coefficient of variation (CV) for the Cq data of all candidate genes. The most stable gene has the lowest SD and CV values. The range of variation in SD should be below 1 (Sarker et al. 2018, Ye et al. 2018).

Results

Primer Specificity and Expression Profiling Analysis of Reference Genes

Gene-specific amplification was confirmed by a melt curve analysis (Supp File 2 [online only]). The melt curves for all genes demonstrated a single peak confirming gene-specific amplification. Agarose gel electrophoresis revealed a single band for all amplified genes (Supp File 1 [online only]). The primers’ regression coefficients (R2) and PCR amplification efficiencies (E) were determined based on the standard curves, showing that R2 values exceeded 0.99, and E ranged from 95 to 116%; the amplification curves were also visibly regular and smooth (Fig. 1). The overall variability of expression for the five candidate reference genes (TUB, 18S, ACT, GAPDH, and EF1α) was studied across nine different temperatures based on an analysis of the raw Cq and SD values (Fig. 2, Table 1). The minimum Cq value of the five reference genes was 15.68 for 18S, indicating it had highest expression abundance, whereas the maximum value was 30.56 for TUB, indicating that it had the lowest expression abundance (Table 1). The CVs for GAPDH, EF1α, TUB, ACT, and 18S were 5.51, 3.25, 9.28, 2.20, and 3.24%, respectively, corresponding to a compact Cq value distribution (Fig. 2), with 18S having the lowest variance in Cq value.

Fig. 1.

Amplification curves and standard curves for the five housekeeping genes.

Fig. 2.

Expression levels of candidate reference genes across different temperature. Distribution of the Cq values obtained for the five housekeeping genes in Calliptamus italicus across nine different temperatures. The boxes indicate the 25th to 75th percentiles. The horizontal lines dividing the boxes indicate the median values, while whiskers represent the maximum and the minimum values. The “+” symbols represent the means.

Expression Stability of Candidate Reference Genes

To evaluate the expression stability of reference genes, the geNorm, BestKeeper, and NormFinder algorithms were used, and the results are described as follows.

geNorm Analysis

The expression stability (M) of each candidate reference gene across the different tested RNA samples was calculated using geNorm (Vandesompele et al. 2002). The genes with the lowest M-values have the most stable expression, and stepwise exclusion of the genes with the highest M-values can, thus, be used to identify combinations of genes with the highest stability. Otherwise, the geNorm algorithm also calculates the average pairwise variation V (V) to determine the optimal number of genes for normalization of RT-qPCR results. When V < 0.15, at least the n + 1 best reference genes are necessary for normalization. In this study, EF1α was predicted to be the most stable gene by geNorm, based on its lowest M-value (1.333), whereas ACT and GAPDH were found to be least stable, with M-values of 1.457 and 1.475, respectively. However, 18S and TUB had variable expression stabilities across the nine different temperature as demonstrated by their M-values >1.5 (Fig. 3A).

Fig. 3.

Evaluation of stability of candidate reference gene expression and pairwise variation V as calculated under different experimental temperature conditions using geNorm software. (A) The gene with the lowest geNorm M-value (<1.5) is considered to have the most stable expression. (B) Pairwise variation (V/V) was calculated between the normalization factors NF and NF to determine the optimal number of reference genes for accurate normalization.

For different temperatures, the recommended optimal numbers of reference genes are shown in Fig. 3B. The V2/3 and V3/4 values were 0.076 and 0.062, respectively, suggesting that the normalization factor should preferably contain three reference genes (Fig. 3B).

NormFinder Analysis

NormFinder analysis assessed the stability value (SV) of each candidate reference gene, with the lowest SV indicating the highest gene expression stability. In this research, the stability ranking results under temperature stress were similar to those of the geNorm analysis. EF1α was the most stable gene, and TUB was the least stable gene across temperature treatments (Fig. 4). Similarly, expression of ACT was more stable than that of GAPDH, in accord with the geNorm analysis results.

Fig. 4.

Evaluation of stability of candidate reference gene expression under different experimental temperature conditions using NormFinder software. The gene with the lowest NormFinder stability value is considered to have the most stable expression.

BestKeeper Analysis

The BestKeeper analysis results summarized in Table 2 indicate that the most stable reference genes were ACT (CV ± SD, 3.12 ± 0.70) and EF1α (CV ± SD, 3.71 ± 0.76) for different temperatures.

Table 2.

Descriptive statistics of five candidate reference genes across the nine different temperatures based on BestKeeper analysis results

Parameters	Gene
	18S	GAPDH	TUB	ACT	EF1
Geometric mean [CP]	15.67	17.82	30.45	22.26	20.59
Arithmetic mean [CP]	15.68	17.83	30.57	22.28	20.60
Min [CP]	14.26	16.17	27.17	20.68	19.31
Max [CP]	16.92	20.39	36.30	24.41	22.43
SD [± CP]	0.61	0.80	2.27	0.70	0.76
CV [% CP]	3.91	4.49	7.43	3.12	3.71

Discussion

Quantitative gene expression analysis across different experimental conditions is necessary for the proper functional analysis of target genes. Further, RT-qPCR is a powerful method to measure gene expression (Pfaffl 2001). However, the accuracy of RT-qPCR in different experimental conditions or species depends strongly on the stability of the reference genes used (Hildyard et al. 2019, Zhao et al. 2019). The ideal reference gene should therefore be stably expressed across a variety of experimental conditions. However, an optimally stable reference gene that meets this criterion across all conditions is almost nonexistent. For example, previous studies across the family Acrididae have shown considerable instability across references genes. The three classical housekeeping genes 18S, GAPDH, and β-actin were influenced by hypobaric hypoxia (Zhao et al. 2012). For adult desert locust brain tissue, the most preferred reference genes were GAPDH, Ubi, and EF1α, but for fifth-instar nymph brain tissue, RP49, EF1α, and ACT have been preferentially recommended as reference genes (Van Hiel et al. 2009). Thus, preliminary experiments are necessary to find appropriate reference genes or gene combinations with stable expression across all experimental conditions (Derveaux et al. 2010, Ye et al. 2018).

In this study, five commonly used reference genes (18S, TUB, EF1α, GAPDH, and ACT) were selected as candidate reference genes for analyses across nine different temperatures. Because this selection of an algorithm for the calculation of the stability of reference genes is inadequate, different analytical approaches based on the geNorm, NormFinder, and BestKeeper algorithms were used.

EF1α is a key protein involved in the elongation cycle of protein biosynthesis, and it is also commonly used as an internal control for normalization for different experiments in insects (Ma et al. 2016, Xu et al. 2017). Similarly, our geNorm and NormFinder algorithm results also show that EF1α was the most stable of the candidate reference genes. Otherwise, based on the geNorm algorithm evaluation, ACT was recommended as the second most stable reference gene across different temperature conditions, in agreement with the NormFinder algorithm. However, there was little disagreement with the BestKeeper algorithm results, as EF1α was chosen as the second ranked reference gene and ACT as the most stable gene. Minimal changes in target gene expression can be masked by the overrepresentation of a reference gene (Raaijmakers et al. 2002, Yang et al. 2018). Although 18S was chosen as the third most stable reference gene, it has a small Cq value (15.68 ± 0.51, Fig. 1, Table 1); therefore, it is not an optimal reference gene. Notably, the geNorm algorithm was employed to evaluate the optimal number of genes for normalization across different temperatures. Our results showed that all reference genes expressed mean pairwise variation values of V2/3 and V3/4 that were <0.15, which indicated that an optimal number of reference genes for normalization in this experiment condition would, therefore, be three. However, although three is the optimal number for calibration in the normalization process, which would make the gene expression profiling more accurate, this is not an absolute standard (Vandesompele et al. 2002). TUB is a common reference gene, and it showed the best stability under various thermal stress conditions in Bemisia tabaci (Dai et al. 2017), but in our study, TUB was the least stable according to geNorm, NormFinder, and BestKeeper algorithms in C. italicus.

This study is the first to systematically analyze reference genes across different temperatures in this species. The most appropriate reference gene in C. italicus is EF1α, but GAPDH and ACT can also be effectively utilized if two normalization genes are needed. These results will enable more accurate quantification of gene expression levels in C. italicus across different temperatures.

Click here for additional data file.

Click here for additional data file.

Click here for additional data file.

Click here for additional data file.

Click here for additional data file.

Click here for additional data file.