JARB Journal of Animal Reproduction and Biotehnology

OPEN ACCESS pISSN: 2671-4639
eISSN: 2671-4663

Article Search

Original Article

Article Original Article

Journal of Animal Reproduction and Biotechnology 2019; 34(4): 280-291

Published online December 31, 2019

https://doi.org/10.12750/JARB.34.4.280

Copyright © The Korean Society of Animal Reproduction and Biotechnology.

Assessment of Suitable Reference Genes for RT-qPCR Normalization with Developmental Samples in Pacific Abalone Haliotis discus hannai

Sang Yoon Lee1, Choul-Ji Park2 and Yoon Kwon Nam1,*

1Department of Marine Bio-Materials and Aquaculture, Pukyong National University, Busan 48513, Korea
2Genetics and Breeding Research Center, National Institute of Fisheries Science, Geoje 53334, Korea

Correspondence to: Yoon Kwon Nam
E-mail: yoonknam@pknu.ac.kr
https://orcid.org/0000-0001-8870-2098

Received: September 28, 2019; Revised: October 25, 2019; Accepted: October 28, 2019

Potential utility of 14 candidate housekeeping genes as normalization reference for RT-qPCR analysis with developmental samples (fertilized eggs to late veliger larvae) in Pacific abalone Haliotis discus hannai was evaluated using four different statistical algorithms (geNorm, NormFinder, BestKeeper and comparative ΔCT method). Different algorithms identified different genes as the best candidates, and geometric mean-based final ranking from the most to the least stable expression was as follow: RPL5, RPL4, RPS18, RPL8, RPL7, UBE2, RPL7A, GAPDH, RPL36, PPIB, EF1A, ACTB and B-TU. The findings were further validated via relative quantification of metallothionein (MT) transcripts using the stable and unstable reference genes, and expression levels of MT were greatly influenced according to the choice of reference genes. In overall, our data suggest that RPL5 and RPS18, either singly or in combination, are appropriate for normalizing gene expression in developmental samples of this abalone species, whereas ACTB, B-TU and EF1A are less stable and not recommended. In addition, our findings propose that standard deviations in geometric ranking as well as geometric mean itself should also be taken into account for the final selection of reference gene(s). This study could be a useful basis to facilitate the generation of accurate and reliable RT-qPCR data with developmental samples in this abalone species.

Keywords: abalone Haliotis discus hannai, developmental samples, reference genes, RT-qPCR assay

Reverse transcription-quantitative polymerase chain reaction (RT-qPCR) is a versatile and routine tool for gene expression analysis (Udvardi et al., 2008; Bustin et al., 2009). However, the accuracy of RT-qPCR results may be affected by both biological and technical variations, and an effective normalization using appropriate reference gene(s) is essential for reliable interpretation of the expression of target genes (Guénin et al., 2009). Although several traditional housekeeping genes such as cytoskeletal β-actin (ACTB) and glyceraldehyde 3’-phosphate dehydrogenase (GAPDH) have been widely used for RT-qPCR assays, a considerable number of reports has also claimed that those traditional reference genes would not always show a stable expression under certain conditions (Li and Shen, 2013; Taylor et al., 2013; Yuan et al., 2014; Lee and Nam, 2016a). Ideally, reference genes must be invariantly expressed and non-regulated by experimental or biological conditions. However, it is widely agreed that no given single reference gene can be universally applied to all experimental conditions and the utility of a reference gene could be largely variable under different conditions (Guénin et al., 2009; Nakayama et al., 2018). Hence, in most situations, it is often unavoidable that the most appropriate reference gene(s) may be empirically determined in a case-specific manner.

Developmental samples (embryos and early larvae) of animals have been often reported to display dynamic modulations of various gene sets related with morphogenesis and organ development (Deschamps and Duboule, 2017; Praggastis and Thummel, 2017). Zygotes of many aquatic animals representing external fertilization should undergo a ‘maternal-zygotic transition’ of gene expression, which makes the interpretation of gene expression profile complex and complicated (Blaxter, 2013; Liu et al., 2014; Romney and Podrabsky, 2017). Newly hatched larvae would be protected from embryonic membrane no longer and they should mount significant changes of gene expression to enter new environment. During early ontogeny, larvae often undergo a process involving conspicuous and relatively abrupt change in body structure, which is called metamorphosis. Metamorphosis is known to be accompanied with expression changes of diverse genes related with cell growth and differentiation (Das et al., 2006; Alves et al., 2016; Zhao et al., 2016). Taken together, careful efforts are needed to precisely interpret the expression pattern of target genes with developmental samples, and undoubtedly the use of appropriate reference genes is a critical component for the accurate quantification of the target gene of interest. Selection of unstable housekeeping genes that are involved in developmental and ontogenic processes may lead to serious misinterpretations.

Pacific abalone, Haliotis discus hannai, is one of the most commercially important mollusk species in Korean aquaculture (Park and Kim, 2013). Due to economic interest, various genetic breeding programs are in progress, including selective breeding, chromosome-set manipulation and interspecific hybridization. Comprehensive understanding genes and its expression involved in the development and ontogeny would be a fundamental requirement for all these breeding investigations with regard to evaluate developmental characteristics and early performances of newly developed breeds. Evaluation of reference genes for RT-qPCR normalization in Haliotis species has been reported with respect to type of tissues and experimental challenges using toxicants and bacteria (Wan et al., 2011; Qiu et al., 2013; López-Landavery, 2014; Lee and Nam, 2016a). However, despite its importance, reference genes for developmental samples of abalone species have not been extensively studied.

The objective of this study was to assess suitable reference genes for RT-qPCR normalization with developmental and larval samples in Pacific abalone, H. discus hannai. In this study, we examined the expression patterns of the 14 housekeeping gene candidates with 8 developmental stages (fertilization to late veliger stages) based on 4 different statistical algorithms (geNorm, NormFinder, BestKeeper and comparative ΔCT method). Further, the influence of reference gene choice on the quantification of a target gene was validated using the stable and unstable reference genes recommended by statistical algorithms.

Abalone specimens and biological sampling

Mature female (n = 5-8) and male (n = 5) abalones were induced to release eggs and sperm based on the conventional method of air exposure followed by ultraviolet (UV)-irradiated seawater treatments. Eggs were washed three times with 1 μm-filtered seawater at 18-19°C and inseminated with sperm using wet methods. Fertilized eggs were placed on static incubator containing 1 μm-filtered seawater at 20°C until hatch. Fertilization rate was estimated with at least 330 randomly chosen embryos as percentage of embryos showing successful progress initial cleavages (at 2 hours post insemination: HPI) out of initial number of eggs inseminated. Hatching success was also calculated as percentage of hatched larvae out of initial number of eggs with 200-300 eggs. Both fertilization rates and hatching success were estimated in triplicates per egg batch. After hatch, about 300,000-500,000 hatchlings were incubated in rectangular tank containing 15 tons of 1 μm-filtered at 19 ± 1°C with a daily water exchange rate of 120%. Rearing of swimming larvae were continued until late veliger stage. Approximately 10,000 embryos were sampled at just fertilized (0 HPI), 2-4 cell stage (2 HPI), 8-16 cell stage (4 HPI), morula stage (5 HPI), gastrula stage (7.5 HPI) and hatching (hatched trocophore; 15 HPI). In addition, about 10,000 swimming larvae at early veliger stage (20 HPH) and late veliger stage (45 HPH) were netted for sampling (Lee and Nam, 2016b). Three independent spawning trials were prepared using different abalone broods in order to prepare three biological replications for each developmental stage. Within each biological replication, triplicate samplings were prepared at each developmental stage as triplicate technical replications.

Candidate reference genes

Candidate reference genes were selected based on previous reports on abalone housekeeping genes (Wan et al., 2011; López-Landavery, 2014; Lee and Nam, 2016a). With cDNA templates from abalone developmental samples, specific amplification of each gene segment was confirmed using end-point RT-PCR and ethidium bromide staining of the RT-PCR band (data not shown). Amplified RT-PCR product were sequenced and correct annotation of each gene was validated with BLAST search against NCBI GenBank (http://www.ncbi.nlm.nih.gov). A total of 14 candidate reference genes were selected: cytoskeletal β-actin (ACTB), β-tubulin (B-TU), elongation factor 1α (EF-1A), glyceraldehyde 3-phosphate dehydrogenase (GAPDH), peptidyl-prolyl cis-trans isomerase B (PPIB), ribosomal protein isoforms L3 (RPL3), L4 (RPL4), L5 (RPL5), L7 (RPL7), L7A (RPL7A), L8 (RPL8), L36 (RPL36), S18 (RPS18), and ubiquitin-conjugating enzyme E2 (UBE2). For GenBank accession numbers, refer to Lee and Nam (2016a). qPCR primers used in this study have 20 bp in length, 50-60% of G + C content, 58-63°C of melting temperature, 156-168 bp in amplicon size.

Nucleic acid preparation and qPCR

Upon being sampled, embryos and larvae were immediately frozen on dry ice and stored at -85°C until used. Total RNA was extracted using Trizol Reagent (Invitrogen; Thermo Scientific, Wilmington, DE, USA) and further purified using RNeasy Plus Micro Kit (Qiagen, Hilden, Germany) according to the manufacturers’ protocols. The quality and quantity of the purified RNA were examined by measuring the absorbance ratios at 260/230 nm and 260/280 nm using a Nano Drop ND1000 spectrophotometer (Thermo Scientific, Wilmington, DE, USA).

A sample of total RNA was reverse transcribed into cDNA using Omniscript RT Kit (Qiagen) according to the manufacturer’s instructions. Synthesized cDNA was 8-fold diluted with sterile distilled water and 2 μL of the diluted cDNA was used as a template for qPCR amplification. qPCR amplification reactions were performed using Light Cycler 480 (Roche Applied Science, Penzberg, Germany) and 2X SYBR Mater Mix (Roche) under the following thermal cycling conditions: 45 cycles of denaturation at 95°C for 5 s, 58°C for 10 s and 72°C for 15 s with an initial denaturation step at 95°C for 2 min. At the end of each run, melting curve analysis was performed from 65°C to 95°C to confirm the specificity of amplification reaction. To validate PCR efficiency (E), standard curve for each gene was prepared with 4-log serial dilution points of a pooled cDNA sample (i.e., a mixture of development and larval cDNAs). The slope, efficiency, and coefficient of determination (R2) were determined with the standard curve using qbaseplus software (Biogazelle, Ghent, Belgium).

Because triplicate samplings were made at each of eight developmental stages from every biological replication, a total of 72 cDNA templates were finally prepared (24 cDNA samples per biological replication). The three samples at each stage from each spawning batch were treated as technical replications. Consequently, within each spawning batch, 336 amplification reactions were carried out to generate 112 median Cq values (8 stage-specific Cq values for each of 14 reference genes). Within this context, a total of 1,008 amplification reactions were made to create 336 Cq values for analysis of expression stability in this study.

Statistical analysis using different software programs

The stability of gene expression was assessed using geNorm (Vandesompele et al., 2002; Hellemans et al., 2007), BestKeeper (Pfaffl et al., 2004), NormFinder (Andersen et al., 2004), and comparative ΔCT method (Silver et al., 2006). Then, an overall final ranking of expression stability for each reference gene was determined by calculating geometric mean value of rankings from the four different algorithms. The geNorm analysis was performed using qbaseplus (Biogazelle). Raw Cq values were converted into relative quantity values (RQ) via the formula RQ = E-ΔCq, where E is the validated amplification efficiency of each gene, and ΔCq is the difference between the Cq value and the minimum Cq of each gene among the samples. Expression stability (geNorm M value) of each candidate reference gene was estimated based on the average pairwise variation for the reference gene with all of the other genes. A low M value represented stable gene expression with the cut-off value of 1.5. The pairwise variation (Vn/n + 1) between the two sequential normalization factors (NFn and NFn + 1) was also estimated to define the optimal number of reference genes needed for developmental and larval samples. geNorm results were compared with three different software algorithms. Microsoft Excel-based BestKeeper calculated stability of reference genes based on the standard deviation (SD) of untransformed Cq values, where the candidate reference gene with the lowest SD is considered as the most stable gene. The NormFinder add-in for Microsoft Excel was used with the same RQ values for geNorm in order to give a stability value with the lowest value as most stable expression. The ΔCT method compares the relative expression of pairs of genes within each sample and calculated the average SD against the other reference genes to confidently identify useful reference genes.

Validation of expression levels based on choice of different reference genes

Effects of choice of different reference genes on relative expression levels of target gene were tested. Based on our previous study, metallothionein (MT) gene was selected as a target gene because it has been reported to represent a significant difference in its mRNA levels between early embryos and swimming larvae (Lee and Nam, 2016b). A segment of MT transcripts (amplicon = 189 bp) was amplified with a primer pair (forward: 5´-GGTACCGACTGCAAGTGTAA-3´ and reverse: 5´-TCATCGGAAGTCATGTGAGC-3´) (Lee and Nam, 2016b). The cDNA templates were chosen from two developmental stages, just fertilized and late veliger. Expression level of MT transcripts in late veliger larvae relative to that in fertilized embryos was normalized with the three most stable and three most unstable reference genes determined in the overall stability ranking. Relative expression levels of MT genes determined by different references were tested using ANOVA (followed by Duncan’s multiple range tests) with a significance level at p < 0.05 and p < 0.01.

Fertilization and hatching success

Fertilization rates of all of three spawning batches were higher than 88% and averaged to be 92.0 ± 2.4% [mean ± SD] (p > 0.05). Hatching success was slightly more variable among batches ranging to 66% to 81% [75.5 ± 6.4%] (p > 0.05) (Fig. 1). No notable sign for the outbreak of mass mortality was observed during the rearing of swimming larvae, although we didn’t quantitatively estimate the larval viability.

Figure 1.

Mean fertilization rate and hatching success of artificially propagated Pacific abalone Haliotis discus hannai batches. Based on ANOVA, statistically significant differences in neither fertilization rate nor hatching success were found among experimental batches used (p > 0.05).



Cq profiles of reference genes

The PCR efficiency of each gene was at least higher than 0.94. The slopes of the standard curves ranged from -3.213 to -3.098, and the correlation coefficients (R2) from 0.991 to 0.998 (Table 1). The melting curve analysis resulted single peak for each gene and no signal was detected in the negative blank (not shown). The quantification cycle (Cq) values were variable among candidate reference genes. In all the replicates, the lowest and highest Cq values were found in ACTB (16.96) and B-TU (29.54). Based on the Cq values, genes showing high expression levels were RPS18, RPL5, EF1A, ACTB, and RPL8 with mean Cq values of 19.99, 20.05, 20.13, 20.39 and 20.41. Conversely, low-abundance genes were RPL36, B-TU, GAPDH, RPL4 and PPIB with mean Cq values of 27.27, 25.87, 25.24, 23.47 and 23,45, respectively (Fig. 2). As shown in Fig. 2, B-TU, ACTB and EF1A showed a relatively wider range of Cq values [% coefficient variations (CV) ranged 6.4 to 11.96], while narrower in RPL7, RPL8 and GAPDH (% CV = 1.01-1.60) (Fig. 2).

Table 1 . Candidate reference genes and qPCR primers used in this study.

Gene symbolFull namePCR efficiency (%)Coefficient of determination (R2) Sequence (5-3)
ACTBCytoskeletal β-actin95.10.991F: GGTATTGTTCTGGACTCTGG
R: GGTGGTGGTGAATGAGTAAC
EF-1AElongation factor 1-alpha95.70.993F: GCTCTCTGGAAGTTTGAGAC
R: CTCCTTCGAGATACCAGCTT
GAPDHGlyceraldehyde-3-phosphate dehydrogenase94.80.991F: ACCGCTACACAGAAGACAGT
R: TACATCAGGTACTGGGACAC
PPIBPeptidyl-prolyl cis-trans isomerase B97.70.997F: CGAGAAAGCAGGACGAATTG
R: AAGTCCCCTCCTTGGATCAT
B-TUTubulin beta102.90.995F: ACATTCACTAGGTGGGGGTA
R: GTACTGACAATGTGGCGTTG
UBE2Ubiquitin-conjugating enzyme E299.50.997F: CCAAGCTCTTCTTAGTGCAC
R: CTCCCCACTTCCATCACTTT
RPL3Ribosomal protein L399.10.991F: TCATTGCACACACCCAGACT
R: CAATGACCTCATCCTGTTCG
RPL4Ribosomal protein L499.70.992F: GCTGCTTCAAGACCGCTTAT
R: TGGCCAGCTTTCTCTGAAAC
RPL5Ribosomal protein L596.80.996F: AGATGAGGATGGCAAACCAG
R: TCGCTGCTCTCAGAGTCAAA
RPL7Ribosomal protein L797.80.994F: CAAGCTGAACACTCCAAACG
R: TCCACAGCACTGATGTTTCC
RPL7ARibosomal protein L7A95.50.995F: GCTGTCGAAAAAGGTTGAGC
R: TGCTTCAAGACAGCGAACTG
RPL8Ribosomal protein L898.50.995F: TGGAAACTACGCCACAGTCA
R: GTCCTGCCTTCAACATTGGT
RPL36Ribosomal protein L36104.10.999F: TGCGCTGAAACTACTTCCTC
R: GACTGTTGTGCTTTCCTCTG
RPS18Ribosomal protein S1896.40.999F: GGTTATACCCGAGAAGTTCC
R: GTCTTTGTGTGCTGACCTCT

Figure 2.

Raw Cq profile of 14 candidate reference genes in all developmental samples in Pacific abalone. In the box plot, the lower and upper boundaries of each box indicate the 25th and 75th percentiles, respectively. A line within the box marks the median value. Whiskers above and below the box indicate the 90th and 10th percentiles, respectively. Outliers are indicated by dots above and/or below the box. Value in parenthesis is the percent coefficient of variation (% CV) for each reference gene.



Time course expression pattern of candidate references during development

Expression modulations in developing embryos and swimming larvae were variable among reference genes. In a broad sense, the patterns of 14 genes were divided into four categories. For Pattern-I, expression levels of reference genes were progressively increased with developments, consequently giving rise to significant differences in expression levels between early embryonic stages and swimming larval stages. Reference genes such as B-TU, ACTB, EF1A and PPIB belonged to this pattern. The Pattern-II could be characterized by gradual decrease of expression levels after early cleavage stages, but after hatch, increase of expression could be observable either at early veliger or late veliger larval stage. This pattern is typically exemplified by RPL7A and GAPDH genes. The Pattern-III is similar, in overall, with the Pattern-II. However, unlike the Pattern-II, rebound or re-increase of expression at swimming larval stages was hardly seen. Typical examples were RPL7 and RPL8. Finally, RPL3 showed a rapid drop of expression level immediately after fertilization that was followed by a quick rebound at morula stage. Afterward, its expression was gradually declined with the progress of development (Pattern-IV) (Fig. 3).

Figure 3.

Expression patterns of 14 candidate reference genes with developmental progress. Time course expressions of reference genes were manually categorized into 4 modulatory patterns based on normalized relative quantity values determined with qbaseplus. Abbreviations of developmental stages are fertilized eggs (FE), early cleavage-1 (EC-1; 2-4 cells), early cleavage-2 (EC-2; 8-16 cells), morula (MO), gastrula (GA), hatched trocophore (HT), early veliger (EV), and late veliger (LV).



geNorm-based analysis of expression stability

With geNorm-based analysis, RPL5 and RPS18 had the lowest M values (0.154 and 0.155, respectively), indicating these two genes are most stably expressed with developmental samples examined in this study. These lowest expression stability values were followed by RPL8 (M = 0.185), RPL7 (0.212) and UBE2 (0.237). Several ribosomal protein genes (RPL36, RPL4, and PRL7A) and GAPDH showed moderate expression stability. Conversely, B-TU (1.041), ACTB (0.847), EF1A (0.596), PPIB (0.500) and RPL3 (0.411) were the five most unstable reference genes in geNorm analysis (Fig. 4A). The optimal number of reference genes determined using geNorm showed a range of V values between 0.038 to 0.155, where the V2/3 value was only 0.081 (Fig. 4B).

Figure 4.

Expression stability of candidate reference genes determined with geNorm algorithm. (A) geNorm M value, (B) geNorm V value. In (B), the proposed cut-off V value (0.15) is indicated with dashed line.



Expression stability based on BestKeeper, NormFinder and comparative ΔCT

BestKeeper identified RPL7 as the most stable gene, and also RPL8, UBE2, and RPL5 as fairly stable references. On the other hand, B-TU, ACTB, EF1A were recognized as the three least stable reference genes. According to NormFinder algorithm, RPL4, RPL7A, PPIB, RPS18 and RPL5 were the five most stable reference genes. Stability ranking of other genes by NormFinder was as follow: RPL3, GAPDH, EF1A, UBE2, RPL8, RPL36, RPL7, ACTB and B-TU. From the comparative ΔCT method, the three most stable genes were RPL4, RPL5 and RPS18, whereas the three most unstable genes were B-TU, ACTB and EF1A (Fig. 5).

Figure 5.

Analysis of expression stability of candidate reference genes determined with BestKeeper, NormFinder and comparative ΔCT method.



Geometric mean-based overall ranking

Expression stability ranking of reference genes determined by individual software programs were shown in Fig. 6A. Different algorithms represented different orders of reference genes according to stability ranking. However, ACTB and B-TU were invariantly identified as the least stable references by all the algorithms (two lowest rankers). Besides these two genes, RPS18, GAPDH and RPL5 showed less variations in rankings across algorithms than did other references, whereas RPL4, PPIB and RPL7 exhibited relatively large fluctuations in ranking depending on algorithms used.

Figure 6.

Stability ranking of 14 candidate reference genes. (A) Ranking determined by individual algorithms (geNorm, NormFinder, BestKeeper and ΔCT). (B) Overall ranking calculated with geometric means of rankings from individual programs. In (B), values in parentheses are geometric means for the most stable (RPL5) and least stable (B-TU), respectively. (C) Calculated geometric standard deviation for the final ranking of each reference gene.



An overall final ranking was calculated by integrating rankings from the four different algorithms (geNorm, BestKeeper, NormFinder, and ΔCT) (Fig. 6B). Based on the geometric mean, most and least recommended references were RPL5 and B-TU, respectively. Ranking from the most to the least stable expression was as follow: RPL5 (geometric mean = 2.51), RPL4 (2.74), RPS18 (3.31), RPL8 (4.16), RPL7 (4.56), UBE2 (5.54), RPL7A (5.58), GAPDH (6.29), RPL36 (7.50), PPIB (7.95), EF1A (10.84), ACTB (13.00) and B-TU (14.00). Calculation of geometric standard deviations for the ranking of each reference showed that ACTB, B-TU and RPS18 received relatively consistent or similar ranks by different algorithms, whereas rankings of RPL7, PPIB and RPL4 were more variable according to algorithms used (Fig. 6C).

Influence of reference gene choice on relative expression of target gene

MT expression of veliger larvae relative to fertilized embryos were measured to be 10.9-, 7.9- and 10.7-fold increase, when normalized with the three most recommended references RPL5, RPL4 and RPS18, respectively. Average fold increase was 9.84 1.26 relative to fertilized embryos. However, normalizations with the three least stable references (EF1A, ACTB and B-TU) showed in significantly different results, in which relative expression level in veliger larvae were 1.6, 0.3 and 0.1 normalized with the reference genes EF1A, ACTB and B-TU, respectively. Calculated average expression level with these three references was 0.70 ± 0.60, indicating 1.4-fold downregulation at veliger larvae relative to fertilized embryos (Fig. 7).

Figure 7.

Relative expression levels of metallothionein transcripts normalized with the three most stable and three least stable reference genes. Statistical differences in expression levels of MT with different normalization references were tested with ANOVA followed by Duncan’s multiple range tests at p < 0.05 and p < 0.01. Means with different letters indicate significant difference.


Biological and technical errors introduced throughout the RT-qPCR process should be normalized using appropriate reference genes by internally controlling for such errors as well as input RNA amounts. High RNA quality and purity are mandatory requirements for reliability and reproducibility of gene expression data (Schmittgen and Livak, 2008; Udvardi et al., 2008). We used only the RNA samples showing the 1.91-2.15 ratios for both 260/280 nm and 260/230. Reasonable PCR efficiencies (higher than 0.94 in this study) with fairly uniform R2 values larger than 0.99 also suggest that quality of RNA and cDNA could be sufficient enough for RT-qPCR assay (Doak and Zair, 2012).

In the present study, we were unable to identify the universal reference gene that were consistently ranked as the best reference across the four software programs. Although the RPL4 was redundantly identified as the most stable gene in NormFinder and ΔCT method, it was only ranked 7th and 8th recommended gene in geNorm (top ranker = RPL5) and BestKeeper (top ranker = RPL7), respectively. This discrepancy may come from different algorithms used in these programs such as ‘pairwise comparison’ or ‘model-based approach’ (Nakayama et al., 2018). To manage this situation, we evaluated overall ranking of each gene based on geometric mean value of individual rankings from the four algorithms (De Santis et al., 2011; Lee and Nam, 2016a; Nakayama et al., 2018). Geometric mean-based ranking identified RPL5 as the best candidate, which was already recommended as the most stable gene by geNorm (lowest M value = 0.154). The RPL5 was ranked 2nd to 5th in other programs. The second gene is RPL4 owing to strong recommendations by NormFinder and ΔCT method as mentioned above. The 3rd recommended gene is the RPS18 with a relative uniformity in rankings across the four different algorithms, although none of program identified it as the best one. Calculation of geometric standard deviation of ranking indicates that the RPS18 (3.64) shows a lower value than all others (5.65-138.78) except ACT and B-TU (both 1.00 corresponding to zero in arithmetic standard deviation) showing invariant ranking irrespective of algorithms. Selection of RPS18 as an internal control for RT-qPCR normalization with developmental samples is congruent with the previous report on other mollusk species (Crassostrea gigas) (Du et al., 2013).

In contrast to variations for top rankers, genes identified to be very unstable represented relatively consistent rankings across software programs. In particular, both B-TU (14th) and ACTB (13th) were always recognized as the two least stable references according to all the four algorithms. Their expression levels are significantly elevated with the progress of development as expressed by normalized RQ values. Fold changes of normalized RQ of B-TU and ACTB between early cleavage and veliger larval stages were more than 15 times, which was apparently dissimilar from other references (maximum 2.3-fold upregulation and 3.4-fold downregulation). Upregulation of these two genes during developmental progress may reflect the increasing demand of cytoskeletons for massive cell growth and morphogenesis particularly in larval period, because swimming abalone larvae undergo dynamic changes of shape and body architecture (Searcy-Bernal et al., 2007; Bunnell et al., 2011; Breuss et al., 2017). Meanwhile, EF1A is an another unstable gene that is consistently ranked 12th by all the algorithms except NormFinder (8th ranked), and changing pattern of normalized RQ during development is similar with those of ACTB and B-TU. GAPDH is one of widely used traditional references (Barber et al., 2005; Zainuddin et al., 2010). However, it was not robustly recommended as a suitable reference gene by any of algorithms (ranked between 7th or 8th in three programs and 4th ranked only in ΔCT method) in this study, which is similar with previous argument on the unsuitability of GAPDH as a reliable reference (Glare et al., 2002; Cho et al., 2008; Sikand et al., 2012). The RPL3 is suggested as the least stable gene out of all ribosomal protein genes tested in this study. Expression pattern of RPL3 during development could be characterized as the downregulation in cleavage stages (2-4 cell stage and further decrease at 8-16 cell stage) relative to fertilized eggs, followed by a quick rebound at morula stage. This change could be explained by a ‘maternal-zygotic transition’ pattern where maternally provisioned transcripts in eggs were replaced with de novo synthesized transcripts by zygotes after early cleavage phase (Harvey et al., 2013; Liu et al., 2014). Such a phenomenon is not clearly observed with other ribosomal protein genes examined in this study, and RPL3-specific fluctuation of expression levels may attribute to its relatively lower rankings than other ribosomal proteins.

Based on the pairwise variation V between two sequential normalization factors containing an increasing number of genes (i.e., geNorm V value), use of two reference genes (V2/3 = 0.081) is believed to be enough for reasonable and accurate normalization in RT-qPCR assay with abalone developmental samples when the proposed cut-off value (V = 0.15) is considered (V = 0.15) (Vandesompele et al., 2002; Hellemans et al., 2007). Finally, we tested the effects of reference choice on the relative expression levels of target gene. The MT gene, which is selected as a target gene model in this study, has been proposed to be closely related with reproduction and gamete quality of marine invertebrates and has also been to be significantly modulated during developmental and ontogenic period in this abalone species (Mao et al., 2012; Lee and Nam, 2016b). Normalization with stable and unstable references resulted in great difference in relative expression levels of MT. Even the normalization with ACTB or B-TU (the two lowest rankers) resulted in the change of MT expression toward opposite direction in contrast to significantly upregulated pattern normalized with highly recommended references RPL5, RPL4 and RPS18. Our finding is undoubtedly suggestive of the importance of selecting appropriate reference genes for the accurate interpretation on the relative expression of target gene in developmental samples (Song et al., 2017). However, despite the overall agreement with choice of top rankers from software programs, we found that normalization of top rankers may not always provide a uniform or similar result regarding the relative expression levels of target gene of interest. In this study, expression level of MT determined with RPL4 (the 2nd gene in the geometric mean ranking) was significantly lower than those with RPL5 (top ranker) and RPS18 (the 3rd ranker). This finding may be associated likely with the large variation in the ranking of RPL4 depending on algorithms used (ranked 1st to 8th positions) in spite of its 2nd top position in the final geometric mean-based ranking. The geometric standard deviation in ranking of RPL4 (43.59) is calculated to be more than 10 times greater than that of the 3rd ranker (RPS18). Hence, this finding suggests that geometric mean-based ranking alone must not be taken as an absolute guide to choose the reference gene, and also that variation or standard deviation in ranking should also be carefully considered to select the final set of reference genes.

In conclusion, potential suitability of 14 candidate housekeeping genes as the normalization reference for RT-qPCR analysis with developmental samples in Pacific abalone was evaluated using four different statistical algorithms. Our data suggest that RPL5 and RPS18, either singly or in combination, are appropriate for normalizing gene expression in developmental samples of this abalone species, whereas ACTB, B-TU and EF1A are less stable and not recommended as references. In addition, our finding also suggests that not only geometric mean-based ranking but also standard deviation in the ranking should be taken into account together for final selection of the best combination of reference gene set. This study could be a useful basis to facilitate the generation of accurate and reliable RT-qPCR data with developmental samples in Pacific abalone.

No potential conflict of interest relevant to this article was reported.

This research was supported by the grant from the Golden Seed Project (GSP), Ministry of Oceans and Fisheries to YKN (#213008-05-2-CG600).

SY Lee, Pukyong Nat’l Univ., Researcher, https://orcid.org/0000-0001-6160-8560

CJ Park, National Institute of Fisheries Science, Researcher, https://orcid.org/0000-0001-8448-6957

YK Nam, Pukyong Nat’l Univ., Professor, https://orcid.org/0000-0001-8870-2098

  1. Alves RN, Gomes SA, Stueber K, Tine M, Thorne MAS, Smaradottir H, Reinhard R, Clark MS, Ronnestad I, and Power DM. (2016) The transcriptome of metamorphosing flatfish. BMC Genomics 17: 413.
    Pubmed KoreaMed CrossRef
  2. Andersen CL, Jensen JL, and Ørntoft TF. (2004) 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. Cancer Res 64: 5245-5250.
    Pubmed CrossRef
  3. Barber RD, Harmer DW, Coleman RA, and Clark BJ. (2005) GAPDH as a housekeeping gene:analysis of GAPDH mRNA expression in a panel of 72 human tissues. Physiol. Genomics 21: 389-395.
    Pubmed CrossRef
  4. Blaxter M. (2013) Development:the maternal-zygotic transition revisited. Curr. Biol 24: R72-75.
    Pubmed CrossRef
  5. Breuss MW, Leca I, Gstrein T, Hansen AH, and Keays DA. (2017) Tubulins and brain development –the origin of functional specification. Mol. Cell. Neuro 84: 58-67.
    Pubmed CrossRef
  6. Bunnell TM, Burbach BJ, Shimizu Y, and Ervasti JM. (2011) β-Actin specifically controls cell growth, migration, and the G-actin pool. Mol. Biol. Cell 22: 4047-4058.
    Pubmed KoreaMed CrossRef
  7. Bustin SA, Benes V, Garson JA, Hellemans J, Huggett J, Kubista M, Mueller R, Nolan T, Pfaffl MW, Shipley GL, Vandesompele J, and Wittwer CT. (2009) The MIQE guidelines:minimum information for publication of quantitative real-time PCR experiments. Clin. Chem 55: 611-622.
    Pubmed CrossRef
  8. Cho YS, Lee SY, Kim KH, and Nam YK. (2008) Differential modulations of two glyceraldehyde 3-phosphate dehydrogenase mRNAs in response to bacterial and viral challenges in a marine teleost Oplegnathus fasciatus (Perciformes). Fish Shellfish Immunol 25: 472-476.
    Pubmed CrossRef
  9. Das B, Cai L, Carter MG, Piao YL, Sharov AA, Ko MSH, and Brown DD. (2006) Gene expression changes at metamorphosis induced by thyroid hormone in Xenopus laevis tadpoles. Dev. Biol 291: 342-355.
    Pubmed CrossRef
  10. De Santis C, Smith-Keune C, and Jerry DR. (2011) Normalizing RT-qPCR data:are we getting the right answers?An appraisal of normalization approaches and internal reference genes from a case study in the finfish Lates calcarifer. Mar Biotechnol 13: 170-180.
    Pubmed CrossRef
  11. Deschamps J, and Duboule D. (2017) Embryonic timing, axial stem cells, chromatin dynamics, and the Hox clock. Genes Dev 31: 1406-1416.
    Pubmed KoreaMed CrossRef
  12. Doak SH, and Zair Z. (2012) Real-time reverse-transcription polymerase chain reaction:technical considerations for gene expression analysis. Genetic Toxicology:Principles and Methods, Methods in Molecular Biology, Parry JM, and Parry EM (eds.) 817, pp.251-270. Springer, NY.
    Pubmed CrossRef
  13. Du Y, Zhang L, Xu F, Huang B, Zhang G, and Li L. (2013) Validation of housekeeping genes as internal controls for studying gene expression during Pacific oyster (Crassostrea gigas) development by quantitative real-time PCR. Fish Shellfish Immunol 34: 939-945.
    Pubmed CrossRef
  14. Glare EM, Divjak M, Bailey MJ, and Walters EH. (2002) β-Actin and GAPDH housekeeping gene expression in asthmatic airways is variable and not suitable for normalising mRNA levels. Thorax 57: 765-770.
    Pubmed KoreaMed CrossRef
  15. Guénin S, Mauriat M, Pelloux J, Van Wuytswinkel O, Bellini C, and Gutierrez L. (2009) Normalization of qRT-PCR data:the necessity of adopting a systematic, experimental conditions-specific, validation of references. J. Exp. Bot 60: 487-493.
    Pubmed CrossRef
  16. Harvey SA, Sealy I, Kettleborough R, Fenyes F, White R, Stemple D, and Smith JC. (2013) Identification of the zebrafish maternal and paternal transcriptomes. Development 140: 2703-2710.
    Pubmed KoreaMed CrossRef
  17. Hellemans J, Mortier G, De Paepe A, Speleman F, and Vandesompele J. (2007) qBase relative quantification framework and software for management and automated analysis of real-time quantitative PCR data. Genome Biol 8: R19.
    Pubmed KoreaMed CrossRef
  18. Lee SY, and Nam YK. (2016a) Evaluation of reference genes for RT-qPCR study in abalone Haliotis discus hannai during heavy metal overload stress. Fish. Aquat. Sci 19: 21.
    CrossRef
  19. Lee SY, and Nam YK. (2016b) Transcriptional responses of metallothionein gene to different stress factors in Pacific abalone (Haliotis discus hannai). Fish Shellfish Immunol 58: 530-541.
    Pubmed CrossRef
  20. Li R, and Shen Y. (2013) An old method facing a new challenge:re-visiting housekeeping proteins as internal reference control for neuroscience research. Life Sci 92: 747-751.
    Pubmed KoreaMed CrossRef
  21. Liu MM, Davey JW, Jackson DJ, Blaxter ML, and Davison A. (2014) A conserved set of maternal genes?Insight from a molluscan transcriptome. Int. J. Dev. Biol 58: 501-511.
    Pubmed KoreaMed CrossRef
  22. Lόpez-Landavery EA, Portillo-Lopez A, Gallardo-Escarate C, and Rio-Portilla MAD. (2014) Selection of reference genes as internal controls for gene expression in tissues of red abalone Haliotis rufescens (Mollusca, Vetigastropoda;Swainson 1822). Gene 549: 258-265.
    Pubmed CrossRef
  23. Mao H, Wang DH, and Yang WX. (2012) Involvement of metallothionein in the development of aquatic invertebrate. Aquatic Toxicol 110-111: 208-2013.
    Pubmed CrossRef
  24. Nakayama T, Okada N, Yoshikawa M, Asaka D, Kuboki A, Kojima H, Tanaka Y, and Haruna S. (2018) Assessment of suitable reference genes for RT-qPCR studies in chronic rhinosinusitis. Sci. Rep 8: 1568.
    Pubmed KoreaMed CrossRef
  25. Park CJ, and Kim SY. (2013) Abalone aquaculture in Korea. J Shellfish Res 32: 17-19.
    CrossRef
  26. Pfaffl MW, Tichopad A, Prgomet C, and Neuvians TP. (2004) Determination of stable housekeeping genes, differentially regulated target genes and sample integrity:Bestkeeper-Excel-based tool using pair-wise correlations. Biotechnol. Lett 26: 509-515.
    Pubmed CrossRef
  27. Praggastis SA, and Thummel CS. (2017) Right time, right place:the temporal regulation of developmental gene expression. Genes Dev 31: 847-848.
    Pubmed KoreaMed CrossRef
  28. Qiu R, Sun B, Fang S, Sun Li, and Liu X. (2013) Identification of normalization factors for quantitative real-time RT-PCR analysis of gene expression in Pacific abalone Haliotis discus hannai. Chin. J. Oceanol. Limnol 31: 421-430.
    CrossRef
  29. Romney AL, and Podrabsky JE. (2017) Transcriptomic analysis of maternally provisioned cues for phenotypic plasticity in the annual killifish Austrofundulus limnaeus. EvoDevo 8: 6.
    Pubmed KoreaMed CrossRef
  30. Schmittgen TD, and Livak KJ. (2008) Analyzing real-time PCR data by the comparative CT method. Nat. Protoc 3: 1101-1108.
    Pubmed CrossRef
  31. Searcy-Bernal R, Perez-Sanchez E, Anguiano-Beltran C, and Flores-Aguilar R. (2007) Metamorphosis and postlarval growth of abalone Haliotis rufescens in a Mexican commercial hatchery. J Shellfish Res 26: 783-787.
    CrossRef
  32. Sikand K, Singh J, Ebron JS, and Shukla GC. (2012) Housekeeping gene selection advisory:glyceraldehyde-3-phosphate dehydrogenase (GAPDH) and β-actin are targets of miR-644a. PLoS One 7: e47510.
    Pubmed KoreaMed CrossRef
  33. Silver N, Best S, Jiang J, and Thein SL. (2006) Selection of housekeeping genes for gene expression studies in human reticulocytes using real-time PCR. BMC Mol. Biol 7: 33.
    Pubmed KoreaMed CrossRef
  34. Song H, Dang X, He YQ, Zhang T, and Wang HY. (2017) Selection of housekeeping genes as internal controls for quantitative RT-PCR analysis of the veined rapa whelk (Rapana venosa). PeerJ 5: e3398.
    Pubmed KoreaMed CrossRef
  35. Taylor DA, Thompson EL, Nair SV, and Raftos DA. (2013) Differential effects of metal contamination on the transcript expression of immune- and stress-response genes in the Sydney Rock oyster, Saccostrea glomerata. Environ Pollt 178: 65-71.
    Pubmed CrossRef
  36. Udvardi MK, Czechowski T, and Scheible WR. (2008) Eleven golden rules of quantitative RT-PCR. Plant Cell 20: 1736-1737.
    Pubmed KoreaMed CrossRef
  37. Vandesompele J, De Preter K, Pattyn F, Poppe B, Van Roy N, De Paepe A, and Speleman F. (2002) Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. Genome Biol 3: research0034.1-0034.11.
    Pubmed KoreaMed CrossRef
  38. Wan Q, Whang I, Choi CY, Lee JS, and Lee J. (2011) Validation of housekeeping genes as internal controls for studying biomarkers of endocrine-disrupting chemicals in disk abalone by real-time PCR. Comp Biochem. Physiol. C 153: 259-268.
    Pubmed CrossRef
  39. Yuan M, Lu Y, Zhu X, Wan H, Shakeel M, Zhan S, Jin BR, and Li J. (2014) Selection and evaluation of potential reference genes for gene expression analysis in the brown planthopper, Nilaparvata lugens (Hemiptera:Delphacidae) using reverse-transcription quantitative PCR. PLoS One 9: e86503.
    Pubmed KoreaMed CrossRef
  40. Zainuddin A, Chua KH, Rahim NA, and Makpol S. (2010) Effect of experimental treatment on GAPDH mRNA expression as a housekeeping gene in human diploid fibroblasts. BMC Mol. Biol 11: 59.
    Pubmed KoreaMed CrossRef
  41. Zhao L, Liu L, Wang S, Wang H, and Jiang J. (2016) Transcriptome profiles of metamorphosis in the ornamented pygmy frog Microhyla fissipes clarify the functions of thyroid hormone receptors in metamorphosis. Sci. Rep 6: 27310.
    Pubmed KoreaMed CrossRef

Article Tools

PDF print Article
Export to Citation Open Access
Google Scholar Send to Email

Share this article on :

Stats or Metrics

Related articles in JARB

more