NORCAL Publications

Research Article

Requirement for Different Normalization Genes for Quantitative Gene Expression Analysis Under Abiotic Stress Conditions in ‘Bangia’ sp. ESS1

Chengze Li ¹, Ryunosuke Irie ¹, Satoshi Shimada ² and Koji Mikami ³*

¹Graduate School of Fisheries Sciences, Hokkaido University, 3-1-1 Minato-cho, Hakodate 041-0821, Japan

²Department of Biology, Ochanomizu University, 2-1-1 Otsuka, Bunkyo-ku, Tokyo 112-8610, Japan

³Faculty of Fisheries Sciences, Hokkaido University, 3-1-1 Minato-cho, Hakodate 041-0821, Japan

*Corresponding author: Koji Mikami, Faculty of Fisheries Sciences, Hokkaido University, 3-1-1 Minato-cho, Japan. Ph: +81-138-40-8899, E-mail: komikami@fish.hokudai.ac.jp

Citation: Mikami K, Li C, Irie R, Shimada S (2019) Requirement for Different Normalization Genes for Quantitative Gene Expression Analysis Under Abiotic Stress Conditions in ‘Bangia’ Sp. ESS1. J Aquat Res Mar Sci 2019: 194-205.

Received Date: 28 June, 2019; Accepted Date: 29 July, 2019; Published Date: 28 August, 2019

Abstract

Quantitative gene expression analysis is indispensable for accurately characterizing the expression levels of genes involved in the response and acclimation of algae to environmental stresses. However, appropriate reference genes in the genus Bangia (order Bangiales) for quantification of stress-inducible gene expression are not well described. Thus, we evaluated the suitability of genes encoding 18S ribosomal RNA (18S rRNA), glyceraldehyde-3-phosphate dehydrogenase (GAPDH), actin, tubulin, elongation factor 1 (EF1), and 60S ribosomal RNA (60S rRNA) as reference genes for real-time quantitative reverse transcription polymerase chain reaction (qRT-PCR) under various stress conditions in ‘Bangia’ sp. ESS1, classified into the ‘Bangia’ 2 group. Our results showed that the 60S rRNA and GAPDH genes were most suitable for gene expression analyses in response to nutrient deficiency, high salinity, and temperature stress, whereas normalization of desiccation-inducible expression required use of the EF1 and tubulin genes. These findings indicate that reference genes for analysis of stress-inducible gene expression by qRT-PCR should be selected according to the environmental stress being studied. Our efforts to identify easily amplified reference genes for various stress responses provide novel insights into the selection of reference genes for seaweeds.

Keywords

Bangia; Gene Expression; Real-Time Quantitative Reverse Transcription Polymerase Chain Reaction (Qrt-PCR); Red Alga; Reference Gene; Stress Response

Introduction

Members of the order Bangiales live in the intertidal zone where organisms periodically experience abiotic stresses such as aberrant temperature, nutrient deficiency, salinity fluctuation, and long periods of desiccation [1-5]. Compared to the recent progress of molecular biological studies on abiotic stress responses in the genus Pyropia [6-11], the study of abiotic stress responses in the genus Bangia is not as advanced. However, given its simple body architecture of gametophytic thalli and ease of a laboratory culture, we recently employed Bangia species for physiological and molecular biological research to elucidate the regulatory mechanisms of the response and acclimation to abiotic stresses in Bangiales. Promotion of the asexual life cycle, a strategy to increase the number of gametophytic clones by producing asexual spores endowed with the gametophyte identity [12], by high-temperature stress was confirmed in Bangia fuscopurpurea [13]. Similar results have been obtained from experiments using B. fuscopurpurea collected from different sites [1,14]. In addition, we have successfully expressed foreign reporter genes transiently in gametophytes of B. fuscopurpurea using a particle bombardment system [15]. Analyses of abiotic-stress-inducible gene expression have also been performed in several laboratories [16-19]. Thus, Bangia species have potential as model red seaweeds for biological research.

Reconstruction of phylogenetic relationships in the genus Bangia has been established recently [20]. For the past few decades, all members of the genus Bangia have been classified into only three species: Bangia atropurpurea, B. fuscopurpurea, and Bangia gloiopeltidicola. However, extensive phylogenetic analysis using sequence information of nuclear small subunit ribosomal RNA (SSU rRNA) and plastid Rubisco large subunit (rbcL) genes from many representatives of Bangia species revealed that the genus Bangia could be divided into four groups: Bangia, ‘Bangia’ 1, ‘Bangia’ 2, and ‘Bangia’ 3. B. atropurpurea and B. gloiopeltidicola fall into groups Bangia and ‘Bangia’ 3, respectively, while B. fuscopurpurea, used in numerous previous studies and thought to represent a single species, is classified into ‘Bangia’ 1 and ‘Bangia’ 2 groups. These findings indicate that specimens previously classified as B. fuscopurpurea are a mixture of phylogenetically close but distinct species. Thus, it is necessary to classify our materials collected at Esashi on the northern Hokkaido Island of Japan, which were described as B. fuscourpurea in previous reports [13,15,21].

Regulation of genes involved in acquisition of stress tolerance is fundamental for intertidal seaweeds to have normal development and growth under strict living conditions [22]. Thus, gene expression analysis is important for understanding the mechanisms of response to environmental stress in intertidal seaweeds. Real-time quantitative reverse transcription polymerase chain reaction (qRT-PCR) allows monitoring of relative changes in gene expression [23] using a normalization strategy that relies on comparison of the target gene with an endogenous control (reference gene). At present, housekeeping genes such as actin, tubulin, and glyceraldehyde-3-phosphate dehydrogenase (GAPDH) as well as ribosomal RNA genes are commonly used as reference genes. These genes are thought to be universally required for basic cellular functions and to be constitutively and stably expressed under different experimental conditions in land plants [24-26]; fungi [27-29]; the green algae Volvox carteri [30], Ulva linza [31], and Chlamydomonas sp. ICE-L [32]; the brown alga Ectocarpus siliculosus [33]; and the red algae Porphyra yezoensis [8,34,35], Pyropia haitanensis [9], and Gracilaria lemaneiformis [36]. Indeed, these housekeeping genes have been employed as reference genes in Chlamydomonas, Volvox, Pyropia, and Porphyra species. However, little is known about appropriate reference genes in Bangia species to date.

In this study, we aimed to identify reference genes suitable for gene expression analysis under various abiotic stresses in Bangia species. Based on our unpublished transcriptome data and information from previous studies [16,18], the six housekeeping genes 18S ribosomal RNA (18S rRNA), GAPDH, EF1, 60S ribosomal RNA (60S rRNA), actin, and tubulin were chosen as candidates for reference genes to evaluate their applicability in quantitative gene expression analysis under four types of stress: desiccation, nutrient deficiency, high salinity, and temperature changes. In addition, a cold-stress-inducible gene encoding delta-12 fatty acid desaturase in B. fuscopurpurea [18] was employed for further validation of the selected reference genes.

Materials and Methods

Algal Strain and Culture Conditions

Thalli of Bangia species used in this study were originally collected at Esashi, Hokkaido, Japan, on May 17, 2010 and maintained clonally in the laboratory as an experimental line. Filamentous gametophytes were cultured in sterilized artificial seawater (SEALIFE; Marinetech, Tokyo, Japan) containing ESS₂ with NaNO₃ as a nitrogen source, vitamins, and trace metal elements [37] under 60µmol photons·m^-2·s^-1 in a short-day photoperiod (10 h light/14 h dark) at 15°C with air filtered through a 0.22- µm filter (Whatman; Maidstone, UK). The culture medium was changed weekly.

Nutrient deficiency, high salinity, desiccation, cold, and heat were selected as abiotic stresses (Table 1). For nutrient deficiency treatment, 0.1g of thalli was treated with ESS₂-less seawater for 24, 48, 72, 96, 120, 144, and 168h. For high-salinity stress, thalli were incubated in medium containing 25, 50 or 100 mM NaCl for 24, 48, and 72h. For desiccation stress, 0.1g of thalli were desiccated in a clean bench until the weight of samples reached 75, 50 or 25% of the initial weight (W_i). To examine recovery from desiccation, thalli weighing 25% of W_i were incubated in normal seawater medium for 1 h. For temperature treatment, 0.1g of thalli were exposed to 5°C or 30°C for 24, 48, and 72 h. All samples were immersed in liquid nitrogen and stored at -80°C prior to RNA extraction.

 

Factor Treatment Duration Desiccation percentage of initial weight 75% 50% 25% 25% seawater 1H Nutrient-deficiency ESS2-less 0D 1D 2D 3D 4D 5D 6D 7D High salinity 25 mM 1d 25 mM 2d 25 mM 3d 50 mM 1d 50 mM 2d 50 mM 3d 100 mM 1d 100 mM 2d 100 mM 3d Temperature changes 5℃ 1d 5℃ 2d 5℃ 3d 30℃ 1d 30℃ 2d 30℃ 3d Table 1: Experimental conditions

 

Phylogenetic Analysis

The rbcL gene of our material was partly amplified with gene-specific primers (5′-AAGTGAACGTTACGAATCTGG-3′ and 5′-GATGCTTTATTTACACCCT-3′) [38] using EX taq polymerase (TaKaRa Bio, Kusatsu, Japan) and sequenced on an ABI Model 3130 Genetic Analyzer (Life Technologies, Carlsbad, USA). The nucleotide sequence of the amplified DNA fragment was deposited inB DDBJ/EMBL/GenBank under accession number MN052802. Specimens and accession numbers of their rbcL gene sequences are listed in (Table 2). A neighbor-joining phylogenetic tree was reconstructed with MEGA 7 software (https://www.megasoftware.net), using ClustalW to align the rbcL nucleotide sequences.

 

Taxon Collection location Accession No. ‘Bangia fuscopurpurea’ BB Bf 1 Bolinas Bay, CA, USA EU289018 ‘Bangia fuscopurpurea’ CMNH UM BF1 Banda, Tateyama, Chiba, Japan HQ687502 ‘Bangia fuscopurpurea’ France Nice, France AF168659 ‘Bangia fuscopurpurea’ Taiwan Taiwan AF168654 ‘Bangia fuscopurpurea’ WA Fisherman's Bay, WA, USA AF169329 ‘Bangia’ maxima Bolinas Bay, CA, USA EU289020 ‘Bangia’ sp. BC Can Ogden Point, Victoria, BC, Canada AF043376 ‘Bangia’ sp. BCH Taylor’s Mistake, Christchurch, South I, NZ HQ687504 ‘Bangia’ sp. BFK Frank Kitts Lagoon, Wellington, North I, NZ HQ687505 ‘Bangia’ sp. BGA Gentle Annie, Westland, South I, NZ HQ687506 ‘Bangia’ sp. BMW Makawhio, North I, NZ HQ687508 ‘Bangia’ sp. BNS Bawley Point, N. of Bateman’s Bay, NSW, Australia HQ687509 ‘Bangia’ sp. BPL Maketu, Bay of Plenty, North I, NZ HQ687510 ‘Bangia’ sp. BRM Kaka Point, Otago, South I, NZ HQ687511 ‘Bangia’ sp. BWP Woodpecker Bay, Paparoa, Wesland, South I, NZ EU570051 ‘Bangia’ sp. CH620 Supseom, Jejudo, Korea HQ728203 ‘Bangia’ sp. MA Woods Hole, MA, USA AF043369 ‘Bangia’ sp. NthBC Can Triple Island, BC, Canada AF043372 ‘Bangia’ sp. OR Lincoln City, OR, USA AF043367 ‘Bangia’ sp. SB Bf 1 Solana Beach, CA, USA EU289019 ‘Bangia’ sp. TX Port Aransas, TX, USA AF043377 ‘Bangia’ vermicularis Golden Gate, San Francisco Bay, CA, USA EU289022 ‘Bangia’ sp. OUCPT-01 Putian, Fujian, China KP747609.1   Total RNA extraction and cDNA synthesis   Table 2: Species and GenBank accession numbers of their rbcL gene sequences used for the     phylogenetic analysis.

 

Total RNA was separately extracted from each sample using a RNeasy Plant Mini kit (Qiagen, Hilden, Germany) and then treated with DNase (TURBO DNA-free TM kit, Invitrogen, Carlsbad, USA) to remove genomic DNA contamination. Purity and concentration of RNA samples were calculated using a GeneQuant pro spectrophotometer (UK), and the integrity of RNA samples was checked using agarose gel electrophoresis. RNA samples of good quality with A₂₆₀/A₂₈₀ ratios ranging from 1.9 to 2.1 were used. First-strand complementary DNA (cDNA) was synthesized from 300 ng of total RNA in a volume of 20µl using a Prime Script 1st strand cDNA synthesis kit (TaKaRa Bio, Kusatsu, Japan) according to the manufacturer’s instructions. The cDNA was diluted 20 times before use as template in qRT-PCR.

Design and Evaluation of Primers for Candidate Reference Genes

Six housekeeping genes encoding 18S rRNA, 60S rRNA, GAPDH, EF1, actin, and tubulin were selected as candidates for reference genes. Primers for qRT-PCR of 18S rRNA and GAPDH genes were available from previous studies [16,18]. Nucleotide sequences of the other four candidate genes were obtained from our unpublished transcriptome analysis: EF1(Unigene64544), actin (CL1634.Contig3), tubulin (CL70.Contig4), and 60S rRNA (Unigene3570). Primers for these genes were designed using Primer Premier 5.0 software (http://www.premierbiosoft.com) with the following conditions: annealing temperature of 60 °C, DNA G+C content of 40 to 60%, and amplicon length of 80 to 200 bp (Table 3). The specificity of primer pairs for each candidate gene to produce PCR fragments of a single expected size was validated by agarose gel electrophoresis of products from PCR reactions performed using Phusion high-fidelity DNA polymerase and GC buffer (Biolabs, Massachusetts, USA) according to the manufacturer’s instructions. A standard curve for each gene was generated using a 5-fold dilution series of one sample over six dilution points measured in three technical replicates, and PCR amplification efficiency was also calculated.

 

 

Quantitative Gene Expression Analysis

Reactions of qRT-PCR were carried out in a total volume of 20µl containing 10µl of 2× SYBR Premix Ex Taq GC, 0.4µl of ROX Reference Dye, 2µl of cDNA template, and 0.4μl (10µM) of each primer, using the SYBR Premix Ex Taq GC kit (TaKaRa Bio, Kusatsu, Japan). The thermal cycling parameters consisted of 95°C for 5 min and 40 cycles of 94 °C for 30 s, 60 °C for 30 s, and 72 °C for 20 s. A dissociation curve was generated by heating from 60 to 95 °C to check for specificity of amplification using an Applied Biosystems 7300 real-time PCR system (Life Technologies, Carlsbad, USA). Each reaction was triplicated, and the data were analysed using one-way ANOVA, where the level of significance was defined at P < 0.05.

 

Evaluation of Gene Expression Stability and Its Validation Using the Delta-12 Fatty Acid Desaturase Gene

Two software packages, geNorm [39] and NormFinder [40], were used to assess the stability of the expression level of each reference gene under different experimental treatment conditions. geNorm ranks the candidate reference genes according to calculation of an expression stability value (M). A lower M value reflects higher stability of gene expression. geNorm also determines the optimal number of reference genes required for accurate normalization of target gene expression by a pairwise variation (V_n/n + 1). Generally, a V_n/n + 1value less than 0.15 indicates no need for an additional reference gene for normalization. NormFinder produces expression stability values (SVs) to evaluate candidate reference genes in a given experimental design and considers intra- and intergroup variations, enabling an estimate of candidate gene stability values for normalization without any influence of coregulated candidate genes. Low SV indicates high stability in normalization.

For validation of these evaluations, the expression of the delta-12 fatty acid desaturase gene was normalized using the reference genes evaluated as suitable by geNorm combined with NormFinder under different abiotic stresses. The primer sequence of delta-12 fatty acid desaturase [18] is indicated in Table S3.

Results

Phylogenetic Classification of The Specimen Used in This Study

A DNA fragment corresponding to the rbcL gene was amplified from our Bangia material, sequenced, and used for phylogenetic analysis. Results indicated that the sample was classified into the ‘Bangia’ 2 group established by Sutherland et al. [20], with the closest relationship to ‘B. fuscopurpurea’ CMNH UM BF1 collected at Banda, Tateyama, Chiba, Japan (Figure 1). We designated our experimental material ‘Bangia’ sp. ESS1 (ESS was derived from Esashi).

 

Figure 1: Phylogenetic identification of the specimen used in the study. Neighbor joining-based phylogenetic tree reconstructed using sequences of rbcL genes from different specimens from the genus Bangia whose DDBJ/EMBL/GenBank accession numbers are listed in Table 2. The specimen used in this study is boxed. Bootstrap values over 50% from 1000 replicates are indicated at nodes. Bar, 0.01 substitutions per site.

 

Amplification Efficiency and Cycle Threshold Values of PCR Reactions Using Candidate Reference Genes

Specificity and efficiency of primers were confirmed by agarose gel electrophoresis of amplicons of the candidate reference genes derived from all cDNA libraries. All the primer pairs amplified single fragments of the expected size (Table 3) that gave rise to one single peak in the melting curve (data not shown). These primer pairs were therefore considered specific and used for evaluation of gene expression levels by qRT-PCR. In addition, the correlation coefficient (R²) and PCR efficiency of each gene were analysed using standard curves. The R² value and efficiency ranged from 0.991 to 0.994 and 91.1% to 105%, respectively (Table 3), both of which were within the acceptable range (0.9910 to 0.9998 for R² and 90–105% for PCR efficiency), indicating the suitability of all primer sets for further gene expression analysis by qRT-PCR.

Cycle threshold (Ct) values represent relative expression levels; low and high values indicate high and low levels of gene expression, respectively. The Ct values of the six candidate reference genes from four different treatments across all samples in ‘Bangia’ sp. ESS1 were calculated as shown in Figure 2. The Ct values of the reference genes varied from 11.5 to 30.7, with the majority between 18 and 28. The 18S rRNA gene showed the highest expression level with Ct values ranging from 11.5 to 13.0, while GAPDH displayed the lowest expression level with Ct values ranging from 26.1 to 30.7. Based on the comparative ranges of Ct values for each reference gene, 18S rRNA showed the most stable gene expression, while 60S rRNA had the most obvious expression variation.

 

Figure 2: Cycle threshold (Ct) values of candidate reference genes across conditions tested. Boxes show the median values (central lines), Q1 (lower outline), and Q3 (upper outline), and whiskers. The whiskers are set at 1.5 times interquartile range (IQR) above Q3 and 1.5 times IQR below Q1. Data were derived from 27 culture samples and repeated for 3 times.

 

However, a simple evaluation of the raw Ct values could not provide sufficient information to determine expression stability of the candidate reference genes. We therefore used geNorm and NormFinder to provide more information on expression stability.

Expression Stability of Candidate Reference Genes

The expression stability values (M) produced by geNorm software allow ranking of the expression stability of each reference gene. Figure 3 ranks the reference genes for effectiveness under each stress treatment based on M, with the least and most stable genes plotted at the left and right, respectively. Other points represent genes with intermediate expression stabilities compared to those of the least and most stable genes. All M values were lower than 1.5, indicating that the expression levels of these genes were relatively stable under all experimental conditions. Among the genes, 60S rRNA and actin were the most stable under nutrient deficiency, salinity, or temperature treatment, whereas EF1 and tubulin showed the best stability under desiccation treatment and conditions combining all treatments (total).

 

Figure 3: Gene expression stability values (M) and ranking of six reference genes in all samples under different treatments calculated using geNorm software. Average expression stability value (M) was calculated by stepwise exclusion of the least stable gene compared across different desiccation, nutrient deficiency, salinity, and temperature conditions, and the combination of all conditions (Total). In each plot, the least stable gene is on the left, and the most stable gene is on the right. Data were derived from 27 culture samples and repeated for 3 times.

 

To provide an independent assessment, the expression stabilities of reference genes were also ranked by NormFinder. Table 4 shows that tubulin and 60S rRNA were the most stable genes under desiccation and temperature treatments, respectively, while GAPDH was highly stable under nutrient deficiency and salinity conditions. GAPDH, tubulin, and EF1 were the most stable genes under conditions combining all treatments.

 

 

The geNorm software also produced information on pairwise variation value (V_n/n+1), which was used to assess the optimal number of reference genes for accurate normalization (Figure 4). In this analysis, a V_n/n+1 value below 0.15 indicates that an additional reference gene is not required. The V_2/3 for each single treatment was below 0.15, implying that use of two stable reference genes is acceptable for normalizing gene expression within each single treatment (Table 5).

 

Figure 4: Determination of the optimal number of reference genes for normalizing gene expression data. Pairwise variation (Vn/n+1), where n represents number of genes) was analysed between normalization factors (NFn andNFn+1) by geNorm software for all different conditions examined (desiccation, nutrient deficiency, salinity, temperature,and all these factors combined). The dotted line at 0.15 is a cut off value, below which an additional reference gene is not needed. Data were derived from 27 culture samples and repeated for 3 times.

 

In addition, V_2/3 and V_3/4 were above 0.15 and below 0.15, respectively, under treatment combining all conditions, revealing that three reference genes are required for gene expression analysis under combined stress conditions (Table 5).

These findings, presented in Figures 3 and 4 and Tables 4 and 5, indicated that use of two or three reference genes ranked as the top two or three in Table 5 is suitable for normalization of stress-inducible gene expression in ‘Bangia’ sp. ESS1. For example, EF1 and tubulin are suitable for desiccation treatment, whereas nutrient deficiency and salinity treatments require use of 60S rRNA, actin, and GAPDH. Combination of all treatments requires EF1, tubulin, and GAPDH as reference genes.

 

 

Validation of reference gene suitability by quantification of expression levels of the delta-12 fatty acid desaturase gene under different stress conditions.

To evaluate the suitability of reference genes selected by our analyses, we quantified the expression levels of the delta-12 fatty acid desaturase gene from ‘Bangia’ sp. ESS1, designated BE1Des12, under normal (15 °C), cold stress (5 °C), and heat stress (30 °C) conditions for 24, 48, and 72 h by qRT-PCR using the most stable (60S rRNA and actin) and least stable (18S rRNA) reference genes. When 60S rRNA, actin, and the combination of both were used as reference genes for normalization, expression levels of BE1Des12 were increased by exposure to cold stress, reaching the highest level at 72 h, while expression levels were not changed under heat stress (Figure 5). By contrast, when the 18S rRNA gene was used as the reference gene, normalized expression levels of BE1Delta12 under both cold and heat stress were not significantly different from those under normal conditions (Figure 5). These experiments confirmed the suitability of actin and 60S rRNA as reference genes for examination of cold-inducible BE1Delta12 expression in ‘Bangia’ sp. ESS1.

 

Figure 5: Gene expression patterns of Delta12-fatty-acid desaturase under temperature stress normalized by the two most stable reference genes, 60S rRNA and actin, and the least stable reference gene, 18S rRNA. Relative expression values shown are mean fold changes compared to samples at 15°C (Control samples, C), with error bars representing standard deviation. Data are from triplicate samples, each with triple technical replicates for qRT-PCR. Different letters above the bars indicate significant differences from the control at 15°C (P < 0.05).

 

These findings also indicated the suitability of the two top-ranked genes in Table 5 as reference genes under the corresponding stress conditions. As BE1Des12 induction by stresses other than cold has not been examined to date, we quantified BE1Des12 expression levels under various stress conditions by normalization with the most suitable reference genes shown in table 5 Desiccation, nutrient deficiency, and high salinity induced expression of BE1Des12, although the kinetics of induction varied among these stresses (Figure 6), suggesting that modifications of membrane fatty acid composition probably occurred following BE1Des12 expression in response to various stresses in ‘Bangia’ sp. ESS1.

 

Figure 6: Gene expression patterns of Delta12-fatty-acid desaturase (BE1Des12) under desiccation, nutrient deficiency, and high salinity, normalized using the two most stable reference genes for each treatment. Relative expression values are mean fold changes compared to control (C) samples, with error bars representing standard deviation. Data are from triplicate samples, each with triple technical replicates for qRT-PCR. Different letters above the bars indicate significant differences from the control at 15°C (P < 0.05). NaCl was used for the high salinity treatment. The label of “25%→1h” in desiccationexperiments means that thalli weighing 25% of the initial weight were incubated in normal seawater medium for 1 h.

 

Discussion

We classified the material used in this study as a member of the ‘Bangia’ 2 group based on phylogenetic analysis (Figure 1). Since this indicated the difficulty of recognizing this material as B. fuscopurpurea, as in our previous study [13], we renamed it ‘Bangia’ sp. ESS1. The plastid genome of a member of the ‘Bangia’ 2 group, designated ‘Bangia’ sp. OUCPT-01 and collected at Putian, Fujian, China, was recently sequenced completely [41]; however, ‘Bangia’ sp. OUCPT-01 is different from ‘Bangia’ sp. ESS because it forms a clade with ‘Bangia’ sp. BGA collected at Gentle Annie, Westland, New Zealand, while ‘Bangia’ sp. ESS1 and ‘B. fuscopurpurea’ CMNH UM BF1 form a distinct sister clade (Figure 1). These findings make the results of biological studies using Bangia species highly confusing because experimental results are partially derived from different Bangia species; thus, it is difficult to obtain a unified view of physiological regulations in a certain species or among species. Therefore, physiological and molecular biological research in Bangia must be performed using several specimens classified phylogenetically and maintained in laboratories to examine biological characteristics of each and among the four ‘Bangia’ groups. This in turn will provide evidence for the current phylogenetic classification of Bangia species.

Our aim here was to identify reference genes suitable for qRT-PCR analysis in studies on the regulation of stress response and tolerance in ‘Bangia’ sp. ESS1. We thus evaluated six candidate reference genes under cold, heat, nutrition deficiency, desiccation, and high salinity, the most common environmental factors influencing growth and development of algae. Stability of candidate gene expression was assessed using two widely used analysis programs, geNorm and NormFinder (Figures 2 and 3; Tables 4 and 5). These programs determined the order of candidate genes in terms of suitability as reference genes under each stress condition and the numbers of reference genes required for accurate normalization. As shown in Tables 4 and 5, GAPDH is a suitable reference gene for gene expression analysis under nutrient deficiency and salinity treatments, whereas normalization of gene expression levels under desiccation stress requires ELF1. GAPDH and ELF1 are necessary for normalization under stress conditions combining all treatments. In addition, actin is a suitable reference for examination of gene expression under nutrient deficiency, temperature changes, and high salinity, while tubulin is a highly stable reference gene under desiccation and combined-stress conditions. Moreover, 60S rRNA is suitable for normalization of gene expression under temperature changes, nutritional deficiency, and high salinity. The suitability of 60S rRNA and actin was indeed validated by quantitative gene expression analysis of low-temperature-inducible expression of BE1Des12 (Figure 5), although Cao et al. [18] employed tubulin and GAPDH genes were employed as references. Therefore, we conclude that there is no perfect reference gene suitable for gene expression analysis under a variety of stress conditions; thus, selection of reference genes appropriate for the desired stress conditions is necessary for performing accurate qRT-PCR analyses in ‘Bangia’ sp. ESS1.

Our findings are in part consistent with previous knowledge. For instance, GAPDH has been established as a suitable reference gene for gene expression studies in the green alga Chlamydomonas sp. ICE-L under different light intensities [42] and in the red alga G. lemaneiformis under temperature stress and at different life cycle stages [36]. EF1 is also recognized as a highly suitable reference gene, showing stable expression, for instance, under the diurnal cycle, high light, high salinity, and UV-B irradiation in the green alga Chlamydomonas sp. ICE-L [42] and under certain stress conditions and at different developmental stages in the red alga P. yezoensis [34]. In addition, the actin gene is useful for studying gene expression responses to temperature changes in G. lemaneiformis [36]; however, it is suitable as a reference gene for studying responses to abiotic stress in P. yezoensis but not gene expression changes at different life cycle stages [34]. The tubulin gene is generally used as a reference gene in the brown alga E. siliculosus, the dinoflagellate Prorocentrum minimum, the green alga U. linza, and the red alga P. yezoensis [31,33,34,43]. Moreover, ribosomal protein genes have been reported as reference genes in the green algae Chlamydomonas sp. ICE-L [32] and V. carteri [44] as well as the diatom Phaeodactylum tricornutum [45]. It is not known how these reference genes were selected in these reports, so we suggest reconfirming the suitability of these reference genes to allow more accurate quantification of stress-inducible gene expression in these algae.

It is notable that 18S rRNA was not suitable as a reference gene for analysis of stress-inducible gene expression in ‘Bangia’ sp. ESS1. The 18S rRNA gene is commonly used as a reference gene in algae [46,47]. For instance, 18S rRNA showed outstanding stability under different light intensities [31] and among life cycle stages [30] in green algae and has been used extensively to normalize gene expression during development and under temperature changes and nutrient deficiency in red algae [9,16]. However, 18S rRNA showed moderate stability across all the treatments we examined in this study. The median Ct value of 12.3 for 18S rRNA was much lower than those for the other reference genes (Figure 2), meaning that the abundance of 18S rRNA mRNA transcripts was much higher. The high abundance of 18S rRNA makes it difficult to reliably subtract the background baseline value when normalizing the mRNA transcripts of target genes. In addition, 18S rRNA did not rank among the three most stable reference genes in any treatment (Table 5). These findings suggest that experiments using 18S rRNA as a reference gene should be revised by identifying more suitable reference genes in algae. For instance, 18S rRNA has been employed as a reference gene for gene expression analyses in the red alga P. yezoensis [10,11,48]; however, it is necessary to find more suitable reference genes for studies of stress-inducible and life cycle generation-dependent gene expression in this organism.

As far as we know, this work is the first to evaluate candidate reference genes for environmental stress-inducible gene expression analysis in Bangia species. Our results indicate that reference genes used in qRT-PCR of stress-inducible gene expression should be employed separately according to differences in environmental stresses to avoid misinterpretation of results in ‘Bangia’ sp. ESS1. The identification of reference genes in this study will facilitate future studies on gene expression in Bangia species to improve our understanding of the molecular mechanisms of responses to environmental stress.

Acknowledgements

Chengze Li was supported by the Ministry of Education, Culture, Sports, Science and Technology of Japan and by the China Scholarship Council.

References

1. Wang WJ, Zhu JY, Xu P, Xu JR, Lin XZ et al. (2008) Characterization of the life history of Bangia fuscopurpurea (Bangiaceae, Rhodophyta) in connection with its cultivation in China. Aquaculture 278: 101-109.
2. Karsten U, West JA (2000) Living in the intertidal zone seasonal effects on heterosides and sunscreen compounds in the red alga Bangia atropurpurea (Bangiales). J Exp Mar Bio Ecol 254: 221-234.
3. Zou D, Gao K (2002) Effects of desiccation and CO2concentrations on emersed photosynthesis in Porphyra haitanensis (Bangiales, Rhodophyta), a species farmed in China. European Journal of Phycology 37: 587-592.
4. Thompson RC, Norton TA, Hawkins SJ (2004) Physical stress and biological control regulate the producer consumer balance in intertidal biofilms. Ecology 85: 1372-1382.
5. Broom JE, Farr TJ, Nelson WA (2004) Phylogeny of the Bangia flora of New Zealand suggests a southern origin for Porphyra and Bangia (Bangiales, Rhodophyta). Mol Phylogenet Evol 31: 1197-207.
6. Gao G, Gao Q, Bao M, Xu J, Li X (2019) Nitrogen availability modulates the effects of ocean acidification on biomass yield and food quality of a marine crop Pyropia yezoensis. Food Chem 271: 623-629.
7. Li C, Kong F, Sun P, Bi G, Li N et al. (2018) Genome wide identification and expression pattern analysis under abiotic stress of mitogen activated protein kinase genes in Pyropia yezoensis. J Appl Phycol 30: 2561-2572.
8. Gao D, Kong F, Sun P, Bi G, Mao Y (2018) Transcriptome wide identification of optimal reference genes for expression analysis of Pyropia yezoensis responses to abiotic stress. BMC Genomics 19: 251.
9. Luo Q, Zhu Z, Zhu Z, Yang R, Qian F et al. (2014) Different responses to heat shock stress revealed heteromorphic adaptation strategy of Pyropia haitanensis (Bangiales, Rhodophyta). PLoS One 9: e94354.
10. Kakinuma M, Nakamoto C, Kishi K, Coury DA, Amano H (2017) Isolation and functional characterization of an ammonium transporter gene, PyAMT1, related to nitrogen assimilation in the marine macroalga Pyropia yezoensis (Rhodophyta). Mar Environ Res 128: 76-87.
11. Uji T, Hirata R, Mikami K, Mizuta H, Saga N (2012) Molecular characterization and expression analysis of sodium pump genes in the marine red alga Porphyra yezoensis. Mol Biol Rep 39: 7973-7980.
12. Takahashi M, Mikami K (2017) Oxidative stress promotes asexual reproduction and apogamy in the red seaweed Pyropia yezoensis. Front Plant Sci 8: 62.
13. Mikami K, Kishimoto I (2018) Temperature promoting the asexual life cycle program in Bangia fuscopurpurea (Bangiales, Rhodophyta) from Esashi in the Hokkaido Island, Japan. Algal Resources 11: 25-32.
14. Notoya M, Iijima N (2003) Life history and sexuality of archeospore and apogamy of Bangia atropurpurea (Roth) Lyngbye (Bangiales, Rhodophyta) from Fukaura and Enoshima, Japan. Fisheries Science 69: 779-805.
15. Hirata R, Takahashi M, Saga N, Mikami K (2011) Transient gene expression system established in Porphyra yezoensis is widely applicable in Bangiophycean algae. Mar Biotechnol (NY) 13: 1038-1047.
16. Wang C, Fan X, Wang G, Niu J, Zhou B (2011) Differential expression of rubisco in sporophytes and gametophytes of some marine macroalgae. PLoS One 6: e16351.
17. Cao TJ, Huang XQ, Qu YY, Zhuang Z, Deng YY, et al. (2017) Cloning and functional characterization of a lycopene betaCyclase from macrophytic red alga Bangia fuscopurpurea. Mar Drugs 15: 116.
18. Cao M, Wang D, Mao Y, Kong F, Bi G et al. (2017) Integrating transcriptomics and metabolomics to characterize the regulation of EPA biosynthesis in response to cold stress in seaweed Bangia fuscopurpurea. PLoS One 12: e0186986.
19. Yokono M, Uchida H, Suzawa Y, Akiomoto S, Murakami A (2012) Stabilization and modulation of the phycobilisome by calcium in the calciphilic freshwater red alga Bangia atropurpurea. Biochim Biophys Acta 1817: 306-311.
20. Sutherland JE, Lindstrom SC, Nelson WA, Brodie J, Lynch MDJ et al. (2011) A new look at an ancient order: Generic revision of the Bangiales (Rhodophyta). Journal of Phycology 47: 1131-1151.
21. Kishimoto I, Ariga I, Itabashi Y, Mikami, K (2019) Heatstress memory is responsible for acquired thermotolerance in Bangia fuscopurpurea. Journal of Phycology.
22. Mikami K (2018) Recent developments in nuclear reverse genetic manipulation that advance seaweed biology in the genomic era. J Aquat Res Mar Sci 1: 01-07.
23. Panina Y, Germond A, Masui S, Watanabe TM (2018) Validation of common housekeeping genes as reference for qPCR gene expression analysis during iPS reprogramming process. Scientific Reports 8: 8716.
24. Czechowski T (2005) Genom-ewide identification and testing of superior reference genes for transcript normalization in Arabidopsis. Plant Physiol 139: 5-17.
25. Tenea GN, Bota PA, Raposo FC, Maquet A (2011) Reference genes for gene expression studies in wheat flag leaves grown under different farming conditions. BMC Res Notes4: 373.
26. Pabuayon IM, Yamamoto N, Trinidad JL, Longkumer T, Raorane ML et al. (2016) Reference genes for accurate gene expression analyses across different tissues, developmental stages and genotypes in rice for drought tolerance. Rice (NY) 9: 32.
27. Kothe E, Cusick KD, Fitzgerald LA, Pirlo RK, Cockrell AL, et al. (2014) Selection and evaluation of reference genes for expression studies with quantitative PCR in the model fungus Neurospora crassa under different environmental conditions in continuous culture. PLoS ONE 9: e112706.
28. Llanos A, Francois JM, Parrou JL (2015) Tracking the best reference genes for RTqPCR data normalization in filamentous fungi. BMC Genomics 16: 71.
29. Chen C, Xie T, Ye S, Jensen AB, Eilenberg J (2016) Selection of reference genes for expression analysis in the entomophthoralean fungus Pandora neoaphidis. Braz J Microbiol 47: 259-265.
30. Kianianmomeni A, Hallmann A (2013) Validation of reference genes for quantitative gene expression studies in Volvox carteri using real-time RT-PCR. Mol Biol Rep 40: 6691-6699.
31. Dong M, Zhang X, Chi X, Mou S, Xu J et al. (2012) The validity of a reference gene is highly dependent on the experimental conditions in green alga Ulva linza. Curr Genet 58: 13-20.
32. Liu C, Wu G, Huang X, Liu S, Cong B (2012) Validation of housekeeping genes for gene expression studies in an ice alga Chlamydomonas during freezing acclimation. Extremophiles 16: 419-425.
33. Le Bail A, Dittami SM, De Franco PO, Rousvoal S, Cock MJ et al. (2008) Normalisation genes for expression analyses in the brown alga model Ectocarpus siliculosus. BMC Mol Biol 9: 75.
34. Kong F, Cao M, Sun P, Liu W, Mao Y (2014) Selection of reference genes for gene expression normalization in Pyropia yezoensis using quantitative real-time PCR. Journal of Applied Phycology 27: 1003-1010.
35. Bhattacharya D, Wu X, Huang A, Xu M, Wang C et al. (2013) Variation of expression levels of seven housekeeping genes at different life history stages in Porphyra yezoensis. PLoS ONE 8: e60740.
36. Ding Y, Sun H, Zhang R, Yang Q, Liu Y, et al. (2014) Selection of reference gene from Gracilaria lemaneiformis under temperature stress. Journal of Applied Phycology 27: 1365-1372.
37. Takahashi M, Saga N, Mikami K (2010) Photosynthesis dependent extracellular Ca²⁺ Influx triggers an asexual reproductive cycle in the marine red macroalga Porphyra yezoensis. American Journal of Plant Sciences 01: 1-11.
38. Ohnishi M, Kikuchi N, Iwasaki T, Kawaguchi R, Shimada S (2013) Population genomic structures of endangered species (CR+EN), Pyropia tenera (Bangiales, Rhodophyta). Japanese Journal of Phycology (Sorui) 61: 87-96.
39. Vandesompele J, De Preter K, Pattyn F, Poppe B, Van Roy, N, et al. (2002) Accurate normalization of real-time quantitative RTPCR data by geometric averaging of multiple internal control genes. Genome Biol 3: RESEARCH0034.
40. Andersen CL, Jensen JL, Orntoft 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.
41. Cao M, Bi G, Mao Y, Li G, Kong F (2018) The first plastid genome of a filamentous taxon ‘Bangia’ sp. OUCPT01 in the Bangiales. Scientific Reports 8: 10688.
42. Mou SL, Zhang XW, Miao JL, Zheng Z, Xu D et al. (2015) Reference genes for gene expression normalization in Chlamydomonas sp. ICE-L by quantitative real-time RT-PCR. Journal of Plant Biochemistry and Biotechnology 24: 276-282.
43. Guo R, Ki JS (2012) Evaluation and validation of internal control genes for studying gene expression in the dinoflagellate Prorocentrum minimum using real-time PCR. Eur J Protistol 48: 199-206.
44. Zhang F, Liu Y, Wang L, Bai P, Ruan L et al. (2018) Molecular cloning and expression analysis of ammonium transporters in tea plants (Camellia sinensis (L.) O. Kuntze) under different nitrogen treatments. Gene 658: 136-145.
45. Siaut M, Heijde M, Mangogna M, Montsant A, Coesel S et al. (2007) Molecular toolbox for studying diatom biology in Phaeodactylum tricornutum. Gene 406: 23-35.
46. Coyne KJ (2010) Nitrate reductase (Nr1) sequence and expression in the harmful alga Heterosigma akashiwo (Raphidophyceae). Journal of Phycology 46: 135-142.
47. Fischer BB, Krieger-Liszkay A, Hideg E, Snyrychova I, Wiesendanger M, et al. (2007) Role of singlet oxygen in chloroplast to nucleus retrograde signaling in Chlamydomonas reinhardtii. FEBS Lett 581: 5555-5560.
48. Kakinuma M, Coury DA, Nakamoto C, Sakaguchi K, Amano H (2008) Molecular analysis of physiological responses to changes in nitrogen in a marine macroalga, Porphyra yezoensis (Rhodophyta). Cell Biol Toxicol 24: 629-639.