https://orcid.org/0000-0002-0663-7366
Figures
Abstract
There has been significant interest in researching the pharmaceutical applications of Industrial hemp since its legalization three years ago. The crop is mostly dioecious and known for its production of phytocannabinoids, flavonoids, and terpenes. Although many scientific reports have showed gene expression analysis of hemp through OMICs approaches, unreliable reference genes for normalization of qRT-PCR data make it difficult to validate the OMICs data. Four software packages: geNorm, NormFinder, BestKeeper, and RefFinder were used to evaluate the differential gene expression patterns of 13 candidate reference genes under osmotic, heavy metal, hormonal, and UV stresses. EF-1α ranked as the most stable reference gene across all stresses, TUB was the most stable under osmotic stress, and TATA was the most stable under both heavy metal stress and hormonal stimuli. The expression patterns of two cannabinoid pathway genes, AAE1 and CBDAS, were used to validate the reliability of the selected reference genes. This work provides useful information for gene expression characterization in hemp and future research in the synthesis, transport, and accumulation of secondary metabolites.
Citation: Deguchi M, Potlakayala S, Spuhler Z, George H, Sheri V, Agili R, et al. (2021) Selection and validation of reference genes for normalization of qRT-PCR data to study the cannabinoid pathway genes in industrial hemp. PLoS ONE 16(12): e0260660. https://doi.org/10.1371/journal.pone.0260660
Editor: Prasanta K. Subudhi, Louisiana State University College of Agriculture, UNITED STATES
Received: July 14, 2021; Accepted: November 13, 2021; Published: December 20, 2021
Copyright: © 2021 Deguchi et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: Our study’s minimal underlying data set was uploaded in Supporting information files.
Funding: RV received grant funding from PA options for Wellness to carry out the research project.
Competing interests: The authors have declared that no competing interests exist.
Introduction
Industrial hemp (Cannabis sativa L.) has a rich history in the civilization of humans because it can provide both phytochemicals and lignocellulosic biomass. This crop originated in Eurasia and is useful all over the world, largely as a fiber crop [1]. Following the emergence of more economically helpful fiber crops, the demand for hemp reduced, and the purpose shifted to usage as a food additive. Hemp seed contains essential fatty acids and proteins and gamma-linolenic acid, which has many health benefits [2]. Hemp seeds and oils are also used to produce nutritional supplements and cosmetics.
Recently, attention has been focused on its rich repertoire of pharmaceutical compounds [3]. Hemp produces a diverse array of phytocannabinoids, terpenes, and phenolic compounds with prominent nutraceutical potential [4]. Among them, phytocannabinoids are the most well-known phytochemicals. The predominant compound, cannabidiol (CBD) and tetrahydrocannabinol (THC) followed by cannabigerol (CBG) and cannabichromene (CBC) are highly promising compounds to improve the quality of human health. They act as therapeutic agents for central nervous system diseases such as epilepsy, inflammation, anxiety, and neurodegenerative disorders such as Parkinson’s, Huntington’s, Tourette’s syndrome, and Alzheimer’s [3]. Terpenes present an array of pharmacological properties, including anxiolytic, antibacterial, anti-inflammatory, and sedative effects on human diseases [5–7].
The increasing popularity of hemp-based phytochemicals has spurred the comprehensive analysis of genome and gene expression studies [8–11]. These analyses are necessary to identify the genes involved in secondary metabolite pathways [12, 13] and have helped to discover transcription factors which are the key proteins that positively or negatively control the synthesis of secondary metabolites [14, 15]. Accurate gene expression studies such as Northern blotting [16], ribonuclease protection assay (RPA) [17], serial analysis of gene expression (SAGE) [18], and quantitative real-time PCR (qRT-PCR) [19] are essential to confirm genomic and transcriptomic data. Among these methods, qRT-PCR is the most frequently used for gene expression analysis because of its high sensitivity, specificity, accuracy, reproductivity, and relatively low cost [20]. qRT-PCR also requires a minimal amount of RNA compared to hybridization-based methods.
The assessment of different samples in the same parameter is evaluated by qRT-PCR [21]. The analysis is used to detect changes in the expression of genes of interest relative to a reference gene. Because of the variances incurred during RNA extraction, DNase treatment, and cDNA synthesis, the reliability of gene expression results can be affected by sample size, RNA degradation, reverse transcription efficiency, and cDNA quality [22]. To provide accurate and reproducible results of gene expression profiles, researchers use reference genes as internal controls. Expression levels vary depending on different environmental conditions, making it critical to identify appropriate and reliable reference genes for each experimental set-up in the respective plant tissue and genotype to prevent biased or misinterpreted data [23–25].
Mangeot-Peter et al. (2014) [26] identified suitable reference genes in hemp stem tissue for accurate expression profiling of cell wall synthesizing genes. Subsequently, Guo et al. (2018) [27] studied seven reference genes in various hemp tissues, such as root, stem, leaf, and flower. To our knowledge, there has been no reports on the suitability of reference genes for normalization of gene expression in hemp under different experimental conditions. The F-box gene is used for hemp gene expression analysis. However, its stability under stress conditions has not been analyzed leading to inaccurate normalization of qRT-PCR analysis [26]. This study aims to evaluate stable reference genes under different abiotic stresses/hormone stimuli in hemp.
Materials and methods
Plant material, greenhouse conditions, generation of clones, growth, and care
Industrial hemp and medical marijuana plants share Cannabis sativa as their common scientific name. Therefore, in this paper, the authors referred to industrial hemp as “hemp”, to distinguish it from medical marijuana.
The hemp strain, Thunderbird, was grown following the approved guidelines for industrial hemp provided by the Pennsylvania Department of Agriculture—Bureau of Plant Industry under the regulated permits IH-16-P-2017 and IH-17-P-2017.
Greenhouse conditions were maintained at 25°C with a 14-hour light photoperiod at 25–40μEm-2s-1. Hemp clones were achieved by collecting a 3-inch segment containing two axillary buds and coating the 45-degree cut with Clonex Rooting gel (Hydrodynamics International, Inc. Lansing, MI). The explant was placed in Root Riot plugs (Hydrodynamics International, Inc. Lansing, MI) and maintained under propagation domes for two weeks at which point they were transferred to four-inch pots containing high porosity soil, HP Mycorrhizae from Pro-Mix (Rivière-du-Loop, Québec, Canada). Genetically identical clones of similar size were obtained by vegetative cuttings from the same female mother plant.
Clones were kept under 24-hour light under propagation domes and 12-hour light during the pre-flowering and flowering periods. The temperature was maintained at 25°C. The humidity for rooting clones was maintained at 65% and decreased gradually to 45% once the clones started to flower. Lost Coast Plant Therapy (Plant Protector, Inc. Loleta, CA) was applied to the clones biweekly at a dilution of 30mL per 4 liters to control pests.
Plant stress treatments
All treatments except UV light treatment were performed in the greenhouse on the same day. Four-week-old, cloned plants grown in small pots were soaked in water including 100mM of mannitol (drought stress), 100mM of NaCl (salt stress), 200μM of CuSO4, 100μM CdCl2, or 100μM of Pb(NO3)2 and 200μM of ZnSO4 (heavy metal stresses), and100μM abscisic acid (ABA), 100μM of methyl jasmonate (MeJA), 1mM of gibberellic acid (GA3), or 100μM of salicylic acid (SA) (hormone treatments) for eight hours. For UV treatment, hemp cloned plants were exposed to UV-C radiation for 10 minutes. After each stress treatment, the 3rd and 4th leaves from the top of the plant were sampled and immediately frozen in liquid nitrogen and stored in the -80°C freezer until total RNA was extracted. All the treatments were performed in three biological replicates. For the mock plants, distilled water was used to soak the hemp plants.
Total RNA extraction and cDNA synthesis
Total RNA was extracted from 100mg of each plant sample using the Spectrum™ Plant Total RNA kit (Sigma Aldrich, St. Louis, MO, USA). RNA concentration and absorbance ratios (A260/280 and A260/230) were measured using a NanoVue Plus spectrophotometer (General Electric Healthcare Limited, UK) to measure the quantity and quality of the total RNA. After treatment with DNase I (TaKaRa Bio, Dalian, China) to remove genomic DNA contamination, 2μg of total RNA was used to synthesize cDNA using the high-capacity cDNA reverse transcription kit (Applied Biosystems, Foster City, CA) according to the manufacturer’s protocol.
Candidate reference genes selection, primer design, and PCR reaction
We identified 13 candidate reference genes (Table 1) and two target genes using a BLAST search from NCBI (https://www.ncbi.nlm.nih.gov/) and the Cannabis sativa genome browser gateway (http://genome.ccbr.utoronto.ca/cgi-bin/hgGateway). The Cannabis sativa genome browser gateway is based on the Purple Kush strain of medical marijuana. Primers were designed based on the sequences of 13 genes using Primer3 Plus (http://www.bioinformatics.nl/cgi-bin/primer3plus/primer3plus.cgi) with the criteria: amplicon size 80–200bp, primer size 18–24bp, Tm 60°C, GC content 45–60%. All primer sequences are listed in Table 1.
qRT-PCR amplification
PCR was performed using cDNA as the template to confirm the specificity of the primers to the target genes. Using a 2% (w/v) agarose gel, all PCR products were analyzed using electrophoresis to confirm a single band of the expected size for each of the primer pairs. To test the PCR amplification efficiency, the regression coefficient (R2) for each gene was calculated using a standard curve generated from a fivefold dilution series of cDNA (1, 1/10, 1/100, 1/1000, and 1/10000) for each primer pair. Based on the slopes of the standard curves, the PCR efficiency of each gene was determined from the respective logarithm of the cDNA dilution and plotted against the mean threshold cycle (Ct) values. The PCR efficiency was calculated using the equation:
where E is the efficiency, and the slope is the gradient of the best-fit line in the linear regression. qRT-PCR was performed with 5μL of SYBR Select Master Mix (Applied Biosystems, Waltham, MA, USA) in a 10μL total reaction mixture containing 400nM of the gene-specific primers and 1μL of cDNA. PCR reaction was performed using a Bio-Rad CFX96 system (Bio-Rad, Hercules, CA, USA) under the following reaction conditions: Initial denaturation at 95°C for 10 minutes, 35 cycles of 95°C for 10 seconds, and 60°C for 1 minute. Three technical replicates were used for each biological replicate and average Ct was used for data analysis. As a negative control, water and total RNA were used instead of cDNA to confirm that there was no amplification from contaminated DNA or hemp genomic DNA.
The stability of reference genes and statistical analysis
Boxplots of quantitative cycle (Cq) values for the 13 candidate reference genes were depicted in all leaf samples with every treatment using the boxplot R package to show the variation of each gene expression. The expression of 13 reference genes was analyzed under 11 different stresses using four algorithms, geNorm, NormFinder, BestKeeper, and RefFinder to rank the stability of the candidate reference genes. The pairwise variation (Vn/Vn+1) between two sequential normalization factors was calculated with geNorm to determine the optimal number of candidate reference genes for accurate normalization [28].
Validation of identified reference genes
Two cannabinoid pathway genes, CBDAS and AAE1, were used as target genes to validate the reliability of the selected reference genes using the most stable candidate reference genes and the least stable reference genes. Primer design and calculation of PCR amplification efficiency for these genes was performed as described above. Relative gene expression levels of CBDAS and AAE1 were calculated using the 2−ΔΔCt method (28). Statistical analysis was performed using a paired t-test (a = 0.05) (28).
Results
PCR specificity and amplification efficiency of the candidate reference genes
Thirteen reference genes (18S, 40S, CHAL, UBE2, EF-1α, F-box, GAD, PCS1, PP2A, SAND, TATA, TIP41, and TUB) were identified from NCBI and the Cannabis sativa genome browser gateway based on a homology search with Arabidopsis genes (Table 1). Primers were designed and used to confirm their specificities based on their amplification efficiency and specificity (Table 1). Single bands were amplified in agarose gel electrophoresis for all the gene primers with predicted sizes (Fig 1). For the qRT-PCR amplification, the PCR efficiency (%) ranged from 91.22 to 113.87, and the regression coefficient (R2) varied from 0.9893 to 0.9994 (Table 1).
Ct values of candidate reference genes
Transcript abundances of 13 candidate reference genes were assessed by qRT-PCR for each gene, tested in triplicates across all 11 treatments and a control, which was 36 biological samples (Fig 2). A majority of the candidate reference genes Ct values ranged from 20 to 30. The lowest expression level with the highest Ct values between 27.4 and 32.1 was PCS1. The EF-1α gene showed the highest expression level, with the lowest Ct values ranging from 18.9 to 25.3. The CHAL gene displayed the highest difference among all 36 samples tested, with a minimum Ct value of 22.3 and a maximum Ct value of 30.9. These Ct value analyses showed that the transcription levels of candidate reference genes are unstable under different stress conditions.
Boxes show the 25th and 75th percentiles, whisker caps represent the minimum and maximum values, lines across the box represent the median Ct-values.
Analysis of reference genes by geNorm
The geNorm was used for evaluating the expression stability of the 13 candidate reference genes (Table 2). Data analysis was calculated based on individual 11 different treatments and three different groups of treatments such as osmotic stress (OS: mannitol, NaCl), heavy metal stress (HM: CdCl2, CuSO4, PbNO3, ZnSO4), and hormonal stimuli (PH: ABA, GA3, MeJA, SA). The total ranking was also shown by combining all 11 treatments together. This algorithm evaluated the stability of reference genes (M) based on the average pairwise variation of all tested genes [29]. In this analysis, the lower the M value, the more stable the gene expression. A reference gene that has an M value less than 1.5 is used for qRT-PCR. PP2A and TIP41 were the most stable reference genes with the lowest M value (0.46) whereas CHAL had an M value of 1.07 and was ranked as the least stable gene. Individually, EF-1a and SAND were the most stably expressed genes under osmotic stresses with an M value of 0.22 while F-box and TATA were the least stably expressed genes. The TUB and TATA genes showed the lowest M values of 0.16 among all of the heavy metal stressed clones. Exposure to hormonal stimuli resulted in PP2A and F-box to be the most stable with an M value of 0.27 and CHAL to be the least stable with an M value of 0.75.
F-box was ranked as the second least stably expressed gene under both osmotic and heavy metal stresses, but it was ranked as the first and second most stable gene under UV and plant hormone treatments, respectively. The TATA gene was least stably expressed under osmotic stresses but was among the top two and three under heavy metal stress and plant hormone stimulus, respectively. CHAL was the least stably expressed in response to UV light application.
geNorm can determine the minimal number of reference genes that should be used to obtain an accurate normalization. The optimal number of reference genes was determined based on the pairwise variation (Vn) between two normalization factors (NFn) composed of an increasing number of reference genes [29]. The threshold value (Vn/Vn+1 = 0.15) indicates if the number of reference genes less than or equal to the value of n is sufficient to use as a reference gene. As shown in Fig 3, the pairwise variation value V2/V3 of all experimental samples was less than 0.15, demonstrating that two reference genes should be sufficient for normalization under all conditions tested.
Analysis of reference genes by NormFinder
NormFinder is a quantity-model-based software and uses complex statistical models to compute the variation between the expression of genes across different biological groups [30]. The lowest expression stability value represents the most stable reference genes. Results from the NormFinder analysis are summarized in Table 3.
The EF-1α and TUB genes were the most stably expressed in all samples and were ranked as fourth and third by geNorm, respectively. The F-box, TATA, and CHAL were ranked as the three least stable genes both by NormFinder and geNorm. The TUB, PCS1, and TATA genes were the most stably expressed under osmotic stress, heavy metal stress, and plant hormone stimuli, respectively. Compared to geNorm, TUB, PCS1, and TATA were ranked as third, sixth, and third positions in each category, respectively. The least stably ex-pressed reference genes under osmotic stress (CHAL, F-box, TATA), heavy metal stress (CHAL, F-box, TATA), and plant hormone stimuli (40S, GAD, CHAL) had similar rankings when compared to geNorm rankings. The GAD and F-box genes were found to be the most stable reference genes under UV stress while PCS1 and CHAL were the least stable. A similar trend was observed in the geNorm analysis.
Analysis of reference genes by BestKeeper
The BestKeeper program is an excel-based algorithm and uses standard deviation (SD) and coefficient of variation (CV) data of the average Ct values for specific treatments [31] (Table 4). Lower CV ± SD values represent higher stability. When using the BestKeeper algorithm, genes with an SD value > 1 are undesirable reference genes [25]. When all samples were taken into consideration, TATA (SD = 0.74), 40S (SD = 0.79), PCS1 (SD = 0.84), EF-1a (SD = 0.90), and TUB (SD = 0.99) were determined to be reliable reference genes. TATA showed the lowest SD among all 13 reference genes in all samples and the SD values were greater than 1 in osmotic stress (1.29) and mannitol (1.82). The 40S gene was ranked as the second most stable candidate in all samples tested, but the SD value of 40S under NaCl and PbNO3 stresses were 1.09 and 1.36, respectively. PCS1 was ranked at the third position in all samples tested and SD values were below 1 in any individual treatment and the three treatment groups. The CHAL gene displayed the highest SD value with 1.95 in all samples displaying that this gene is not ideal for gene expression normalization. The CHAL gene exhibited an SD value less than 1 only under GA3 (0.71), ABA (0.27), CdCl2 (0.58), and CuSO4 (0.98) treatments.
Data after gene symbols mean Std ± CV%.
Analysis of reference genes by RefFinder
RefFinder is a web-based tool for comprehensive analysis that integrates geNorm, NormFinder, Delta Ct, and BestKeeper approaches [32]. The reference genes were ranked from the most stable (M value is the lowest) to the least stable expression (M value is the highest) using RefFinder (Table 5). Among them, the most stable candidate was the EF-1α gene, followed by the TUB gene in all samples. The EF-1α and TUB genes were also ranked in third and first places under osmotic stress conditions, respectively. The TATA gene was most stably expressed under heavy metal and plant hormone treatments while this gene was the least stable under osmotic stress. The CHAL gene was ranked as the least stable gene in all samples tested. The GAD and CHAL genes were the most and least stably expressed genes respectively under UV application, which was the same findings as to the NormFinder software.
Validation of selected reference genes
To validate the selected reference genes, gene expression levels of AAE1 and CBDAS were measured (Fig 4). Each of the two most stable reference genes, EF-1α and TUB, a combination of these two stable reference genes (EF-1α+TUB), and the least stable reference gene (CHAL) were used as internal controls. AAE1 expression was significantly reduced under drought (Mannitol) and salinity (NaCl) stresses. EF-1α, TUB, and a combination of EF-1α and TUB were used for normalization of qRT-PCR analysis. There was no significant difference in the AAE1 expression between the mock treatment and osmotically stressed samples (Mannitol and NaCl) when CHAL was used as an internal control. The expression of CBDAS was also reduced under osmotic stresses when expression data was normalized with EF-1α, TUB, and a combination of EF-1α and TUB unless the CBDAS expression under NaCl stress was normalized with the TUB gene. When CHAL was used as a reference gene, CBDAS gene expression was reduced under mannitol treatment but there was no difference between mock and NaCl treatments.
Error bars for qRT-PCR show the standard error of three replicates for EF-1α, TUB, and CHAL, and six replicates for a combination of EF-1α and TUB. The asterisk represents that there is a significant difference in the comparison with the mock treatment by the statistical analysis (P < 0.05) in paired T-tests.
Discussion
Industrial hemp is from the plant species Cannabis sativa and has gained importance as a medicinal crop because of its potential to produce secondary metabolites such as cannabinoids, terpenes, and phenolic compounds [33]. According to Schluttenhofer and Yuan (2017) [4], hemp was cultivated for commercial or research purposes in at least 47 countries in 2017 and the global hemp market doubled from the year 2016 to 2020. Recently, a comprehensive gene expression analysis is aimed at elucidating the metabolic pathways for cannabinoids and terpene synthesis to improve hemp traits [34–36]. To validate this data, qRT-PCR analysis is suitable, however, appropriate hemp reference genes for accurate gene expression analysis have not been well established.
In this report, we evaluated 13 hemp reference genes under 11 different stress conditions. Research in other plant species has revealed that different environmental conditions would require unique reference genes to accurately interpret expression levels [37, 38]. Eleven different conditions including osmotic stresses, heavy metal stresses, plant hormone stimulus, and UV light application were reported to affect the cannabinoid synthesis [28, 39–41]. The results obtained from geNorm, NormFinder, BestKeeper, and RefFinder were not consistent, particularly BestKeeper which was much more distinct from the other software methods (Tables 2–5). This finding was expected because the BestKeeper algorithm evaluates data differently when compared to the three other programs [42].
To rank the most suitable reference genes across all treatments, there was no unanimity when compared to four different algorithms (Table 6), which represented the combined results obtained from four programs. In most cases, one candidate gene was ranked as the most stable gene by two or three programs, which showed that it might be a good reference gene under various treatments. Based on the combined rankings of the four programs used in our study, the overall results showed that the most stable genes varied while the least stable genes were almost the same. Across all plants tested, both NormFinder and RefFinder determined EF-1α as the most stable gene in all samples tested. In previous reports, EF-1α was demonstrated to be the most stable gene under different stresses in a variety of crops such as tobacco [43], maize [44], soybean [45], potato [46, 47]. Interestingly, this gene was not the most stable in any of the three groups (OS, HM, PH). The TUB gene appeared to be best the candidate under osmotic stresses because this gene was ranked as the most stable by both NormFinder and RefFinder which is consistent with the results obtained in Parsley under abiotic stresses [37]. Under heavy metal stress, TATA was ranked as the most stable gene by BestKeeper and RefFinder and the second most stable gene by geNorm. This gene was also ranked as the best reference gene in hormone stimuli by NormFinder and RefFinder. Interestingly, TATA was the least stable gene under osmotic stresses by geNorm, NormFinder, and RefFinder. TUB was the most stably expressed gene under osmotic stresses, whereas TATA was ranked as the best stable gene under both heavy metal stress and hormone stimuli. Unlike most stable genes, CHAL was found to be the least stable gene in most of the rankings with all samples and the three treatment groups (OS, HM, PH) when analyzed by all four programs. According to Wang et al. 2015 [48], candidate genes showing a high level of variation of Ct values should be avoided as internal controls. Our results showed that variation of the Ct value in CHAL was highest among all 13 reference genes (Fig 2), which is consistent with the fact that CHAL was ranked as the least stable by all four programs used in this study.
In previous Cannabis qRT-PCR studies, the F-box gene has been used as an internal control for qRT-PCR [28, 49, 50]. Mangeot-Peter et al. (2016) [26] performed the reference gene analysis in hemp stems and concluded that the F-box gene was ranked as one of the most stable genes and Histone 3 as the least stable gene among 12 reference genes tested under normal conditions. In this study, however, the F-box gene was the second least stable gene by RefFinder and the third least stable gene determined by both the geNorm and NormFinder programs when all samples were analyzed. Based on our group rankings (OS, HM, PH), F-box was ranked the second least stable genes by geNorm, BestKeeper, and NormFinder under both osmotic and heavy metal stresses. The F-box gene was stably expressed under normal conditions in hemp leaves (S2 Fig) and relatively stable under hormone stimuli as evident by its second position as ranked by both geNorm and RefFinder. These results show that F-box may not be a suitable reference gene for hemp qRT-PCR analysis under osmotic and heavy metal stresses. However, it could be acceptable as a reference gene under normal and plant hormone treatments. Overall, our study suggests that the F-box gene may not be the best reference gene for C. sativa, particularly in plant stress-related studies. Guo et al. (2018) [27] have studied the stability of reference genes in different hemp tissues/organs. They ranked ubiquitin and EF-1α as the most stable genes in leaf samples at different stages, and PP2A as the least stable gene in different organs. Notably, EF-1α was the most stable reference gene in our global ranking, showing that EF-1α is most stable under the normal condition and different abiotic stresses and hormonal stimuli.
Many studies have proved that the use of more than one reference gene enables the possibility of avoiding variations and achieving more accurate normalization of qPCR data [29]. To assess the optimal number of reference genes for the normalization of qRT-PCR data, we used the geNorm program to perform a stepwise calculation of the pairwise variation (Vn/Vn+1) between sequential normalization factors. In this analysis, a Vn/Vn+1 < 0.15 indicates that introducing an additional reference gene for normalization is unnecessary. Under all treatments, V2/V3 values were less than 0.15, which indicated that two reference genes were enough for the normalization of the real-time PCR data under any treatments in this study.
To validate the reliability of the selected reference genes, we measured the relative expression of two cannabinoids pathway genes using EF-1α and TUB as the most stable reference genes and CHAL as the least stable reference gene (Fig 4). Since CBDA content is decreased by the influence of osmotic stress [41], we measured the expression of AAE1 and CBDAS genes that are involved in the rate-determining enzymatic reactions leading to CBDA synthesis under drought and salinity stresses [51, 52]. The expression of these two genes was significantly reduced under drought and salinity stresses when qRT-PCR data were normalized by EF-1α, TUB, and the combination of EF-1α and TUB. Notably, the expression level of both genes was normalized by CHAL under salinity stress and did not show a significant difference when compared with mock plants. These results suggest that EF-1α and TUB genes individually or in combination are suitable reference genes for hemp under osmotic stresses. Our validation study demonstrated the effectiveness of the ranking of reference genes by the programs used geNorm, NormFinder, and RefFinder.
To the best of our knowledge, this study is the first report that performed a systematic analysis of hemp reference genes under different abiotic stresses and hormonal stimuli. The knowledge obtained in this study could contribute to enhancing future hemp research related to the elucidation of mechanisms involved in the synthesis, transport, and accumulation of abundant secondary metabolites in hemp.
Supporting information
S1 Fig. Original gel picture for the amplification of 13 candidate genes using cDNA synthesized from hemp leaf sample to confirm primer specificity and amplicon size.
https://doi.org/10.1371/journal.pone.0260660.s001
(TIF)
S2 Fig. Expression level variability of each candidate reference gene to examine among different tissues (n = 12).
https://doi.org/10.1371/journal.pone.0260660.s002
(TIF)
S3 Fig. Hemp pictures before stress treatments.
https://doi.org/10.1371/journal.pone.0260660.s003
(PDF)
S4 Fig. Hemp pictures after stress treatments.
https://doi.org/10.1371/journal.pone.0260660.s004
(PDF)
S1 Table. Ct values to make a standard curve generated from different dilution series of cDNA for each primer pair.
https://doi.org/10.1371/journal.pone.0260660.s005
(XLSX)
S2 Table. Ct values to study expression level variability of each candidate reference gene under different stresses in all leaf samples.
https://doi.org/10.1371/journal.pone.0260660.s006
(XLSX)
S3 Table. Ct values to analyze the reference genes using four algorithms, geNorm, NormFinder, BestKeeper, and RefFinder.
https://doi.org/10.1371/journal.pone.0260660.s007
(XLSX)
S4 Table. Ct values to evaluate the validation of identified reference genes.
https://doi.org/10.1371/journal.pone.0260660.s008
(XLSX)
S5 Table. Ct values to study expression level variability of each candidate reference gene in different tissues.
https://doi.org/10.1371/journal.pone.0260660.s009
(XLSX)
Acknowledgments
We would like to thank Dr. Yuka Imamura, Ph.D. of Penn State College of Medicine Genome Sciences and Bioinformatics Facility for data analysis services. The authors would like to thank Ms. Jessica Wolf for proofreading the manuscript.
References
- 1. Frassinetti S, Moccia E, Caltavuturo L, Gabriele M, Longo V, Bellani L, et al. Nutraceutical potential of hemp (Cannabis sativa L.) seeds and sprouts. Food Chemistry. 2018;262:56–66. pmid:29751921
- 2. Bonini SA, Premoli M, Tambaro S, Kumar A, Maccarinelli G, Memo M, et al. Cannabis sativa: A comprehensive ethnopharmacological review of a medicinal plant with a long history. Journal of Ethnopharmacology. 2018;227:300–15. pmid:30205181
- 3. Izzo AA, Borrelli F, Capasso R, Di Marzo V, Mechoulam R. Non-psychotropic plant cannabinoids: new therapeutic opportunities from an ancient herb. Trends in Pharmacological Sciences. 2009;30(10):515–27. pmid:19729208
- 4. Schluttenhofer C, Yuan L. Challenges towards Revitalizing Hemp: A Multifaceted Crop. Trends in Plant Science. 2017;22(11):917–29. pmid:28886910
- 5. Chadwick B, Miller M, Hurd Y. Cannabis Use during Adolescent Development: Susceptibility to Psychiatric Illness. 2013;4(129).1–8.
- 6. McPartland JM, Russo EB. Cannabis and Cannabis Extracts. Journal of Cannabis Therapeutics. 2001;1(3–4):103–32.
- 7. Vasas A, Hohmann J. Euphorbia Diterpenes: Isolation, Structure, Biological Activity, and Synthesis (2008–2012). Chemical Reviews. 2014;114(17):8579–612. pmid:25036812
- 8. Braich S, Baillie RC, Jewell LS, Spangenberg GC, Cogan NOI. Generation of a Comprehensive Transcriptome Atlas and Transcriptome Dynamics in Medicinal Cannabis. Scientific Reports. 2019;9(1):16583. pmid:31719627
- 9. McGarvey P, Huang J, McCoy M, Orvis J, Katsir Y, Lotringer N, et al. De novo assembly and annotation of transcriptomes from two cultivars of Cannabis sativa with different cannabinoid profiles. Gene. 2020;762:145026. pmid:32781193
- 10. van Bakel H, Stout JM, Cote AG, Tallon CM, Sharpe AG, Hughes TR, et al. The draft genome and transcriptome of Cannabis sativa. Genome Biology. 2011;12(10):R102. pmid:22014239
- 11. Zager JJ, Lange I, Srividya N, Smith A, Lange BM. Gene Networks Underlying Cannabinoid and Terpenoid Accumulation in Cannabis. Plant Physiol. 2019;180(4):1877–97. pmid:31138625
- 12. Booth JK, Bohlmann J. Terpenes in Cannabis sativa—From plant genome to humans. Plant science: an international journal of experimental plant biology. 2019;284:67–72. pmid:31084880
- 13. Gülck T, Møller BL. Phytocannabinoids: Origins and Biosynthesis. Trends Plant Sci. 2020;25(10):985–1004. pmid:32646718
- 14. Bassolino L, Buti M, Fulvio F, Pennesi A, Mandolino G, Milc J, et al. In Silico Identification of MYB and bHLH Families Reveals Candidate Transcription Factors for Secondary Metabolic Pathways in Cannabis sativa L. Plants (Basel). 2020;9(11):1540. pmid:33187168
- 15. Marks MD, Tian L, Wenger JP, Omburo SN, Soto-Fuentes W, He J, et al. Identification of candidate genes affecting Delta9-tetrahydrocannabinol biosynthesis in Cannabis sativa. J Exp Bot. 2009;60(13):3715–26. pmid:19581347
- 16. Kevil CG, Walsh L, Laroux FS, Kalogeris T, Grisham MB, Alexander JS. An Improved, Rapid Northern Protocol. Biochemical and Biophysical Research Communications. 1997;238(2):277–9. pmid:9299493
- 17. Stalder AK, Pagenstecher A, Kincaid CL, Campbell IL. Analysis of Gene Expression by Multiprobe RNase Protection Assay. Methods in molecular medicine. 1999;22:53–66. pmid:21380823
- 18. Velculescu VE, Zhang L, Vogelstein B, Kinzler KW. Serial analysis of gene expression. Science (New York, NY). 1995;270(5235):484–7. pmid:7570003
- 19. Nolan T, Hands RE, Bustin SA. Quantification of mRNA using real-time RT-PCR. Nature Protocols. 2006;1(3):1559–82. pmid:17406449
- 20. Taylor SC, Nadeau K, Abbasi M, Lachance C, Nguyen M, Fenrich J. The Ultimate qPCR Experiment: Producing Publication Quality, Reproducible Data the First Time. Trends in Biotechnology. 2019;37(7):761–74. pmid:30654913
- 21. Deepak S, Kottapalli K, Rakwal R, Oros G, Rangappa K, Iwahashi H, et al. Real-Time PCR: Revolutionizing Detection and Expression Analysis of Genes. Current genomics. 2007;8(4):234–51. pmid:18645596
- 22. Boulter N, Suarez FG, Schibeci S, Sunderland T, Tolhurst O, Hunter T, et al. A simple, accurate and universal method for quantification of PCR. BMC Biotechnology. 2016;16(1):27. pmid:26956612
- 23. Ambroise V, Legay S, Guerriero G, Hausman JF, Cuypers A, Sergeant K. Selection of Appropriate Reference Genes for Gene Expression Analysis under Abiotic Stresses in Salix viminalis. Int J Mol Sci. 2019;20(17). pmid:31466254
- 24. Niu X, Chen M, Huang X, Chen H, Tao A, Xu J, et al. Reference Gene Selection for qRT-PCR Normalization Analysis in kenaf (Hibiscus cannabinus L.) under Abiotic Stress and Hormonal Stimuli. Front Plant Sci. 2017;8:771. pmid:28553304
- 25. Yu Y, Zhang G, Chen Y, Bai Q, Gao C, Zeng L, et al. Selection of Reference Genes for qPCR Analyses of Gene Expression in Ramie Leaves and Roots across Eleven Abiotic/Biotic Treatments. Scientific Reports. 2019;9(1):20004. pmid:31882847
- 26. Mangeot-Peter L, Legay S, Hausman J-F, Esposito S, Guerriero G. Identification of Reference Genes for RT-qPCR Data Normalization in Cannabis sativa Stem Tissues. International Journal of Molecular Science. 2016;17(9):1556. pmid:27649158
- 27. Guo R, Guo H, Zhang Q, Guo M, Xu Y, Zeng M, et al. Evaluation of reference genes for RT-qPCR analysis in wild and cultivated Cannabis. Bioscience, Biotechnology, and Biochemistry. 2018;82(11):1902–10. pmid:30130459
- 28. Husain R, Weeden H, Bogush D, Deguchi M, Soliman M, Potlakayala S, et al. Enhanced tolerance of industrial hemp (Cannabis sativa L.) plants on abandoned mine land soil leads to overexpression of cannabinoids. PLoS One. 2019;14(8):e0221570. pmid:31465423
- 29. Vandesompele J, De Preter K, Pattyn F, Poppe B, Van Roy N, De Paepe A, et al. Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. Genome Biology. 2002;3(7):research0034.1–0034.11. pmid:12184808
- 30. Andersen CL, Jensen JL, Ørntoft TF. 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 research. 2004;64(15):5245–50. pmid:15289330
- 31. Pfaffl MW, Tichopad A, Prgomet C, Neuvians TP. Determination of stable housekeeping genes, differentially regulated target genes and sample integrity: BestKeeper—Excel-based tool using pair-wise correlations. Biotechnology Letters. 2004;26(6):509–15. pmid:15127793
- 32. Xie F, Xiao P, Chen D, Xu L, Zhang B. miRDeepFinder: a miRNA analysis tool for deep sequencing of plant small RNAs. Plant Molecular Biology. 2012;80(1):75–84. pmid:22290409
- 33. Andre CM, Hausman J-F, Guerriero G. Cannabis sativa: The Plant of the Thousand and One Molecules. Front Plant Sci. 2016;7:19-. pmid:26870049
- 34. Adal AM, Doshi K, Holbrook L, Mahmoud SS. Comparative RNA-Seq analysis reveals genes associated with masculinization in female Cannabis sativa. Planta. 2021;253(1):17. pmid:33392743
- 35. Gao C, Cheng C, Zhao L, Yu Y, Tang Q, Xin P, et al. Genome-Wide Expression Profiles of Hemp (Cannabis sativa L.) in Response to Drought Stress. International Journal of Genomics. 2018;3057272. pmid:29862250
- 36. Prentout D, Razumova O, Rhoné B, Badouin H, Henri H, Feng C, et al. An efficient RNA-seq-based segregation analysis identifies the sex chromosomes of Cannabis sativa. Geome Res. 2020; 30(2):164–172. pmid:32033943
- 37. Li M-Y, Song X, Wang F, Xiong A-S. Suitable Reference Genes for Accurate Gene Expression Analysis in Parsley (Petroselinum crispum) for Abiotic Stresses and Hormone Stimuli. Front Plant Sci. 2016;7:1481. pmid:27746803
- 38. Poli M, Salvi S, Li M, Varotto C. Selection of reference genes suitable for normalization of qPCR data under abiotic stresses in bioenergy crop Arundo donax L. Scientific Reports. 2017;7(1):10719. pmid:28878356
- 39. Burgel L, Hartung J, Schibano D, Graeff-Hönninger S. Impact of Different Phytohormones on Morphology, Yield and Cannabinoid Content of Cannabis sativa L. Plants (Basel). 2020;9(6):725. pmid:32521804
- 40. Eichhorn Bilodeau S, Wu B-S, Rufyikiri A-S, MacPherson S, Lefsrud M. An Update on Plant Photobiology and Implications for Cannabis Production. Front Plant Sci. 2019;10:296. pmid:31001288
- 41. Yep B, Gale NV, Zheng Y. Aquaponic and Hydroponic Solutions Modulate NaCl-Induced Stress in Drug-Type Cannabis sativa L. Front Plant Sci. 2020;11:1169. pmid:32849724
- 42. Guo J, Ling H, Wu Q, Xu L, Que Y. The choice of reference genes for assessing gene expression in sugarcane under salinity and drought stresses. Scientific Reports. 2014;4(1):7042. pmid:25391499
- 43. Schmidt GW, Delaney SK. Stable internal reference genes for normalization of real-time RT-PCR in tobacco (Nicotiana tabacum) during development and abiotic stress. Molecular genetics and genomics: MGG. 2010;283(3):233–41. pmid:20098998
- 44. Lin Y, Zhang C, Lan H, Gao S, Liu H, Liu J, et al. Validation of Potential Reference Genes for qPCR in Maize across Abiotic Stresses, Hormone Treatments, and Tissue Types. PLOS ONE. 2014;9(5):e95445. pmid:24810581
- 45. Saraiva KD, Fernandes de Melo D, Morais VD, Vasconcelos IM, Costa JH. Selection of suitable soybean EF1α genes as internal controls for real-time PCR analyses of tissues during plant development and under stress conditions. Plant cell reports. 2014;33(9):1453–65. pmid:24820128
- 46. Nicot N, Hausman J-F, Hoffmann L, Evers D. Housekeeping gene selection for real-time RT-PCR normalization in potato during biotic and abiotic stress. J Exp Bot. 2005;56(421):2907–14. pmid:16188960
- 47. Tang X, Zhang N, Si H, Calderón-Urrea A. Selection and validation of reference genes for RT-qPCR analysis in potato under abiotic stress. Plant methods. 2017;13:85. pmid:29075311
- 48. Wang H-L, Li L, Tang S, Yuan C, Tian Q, Su Y, et al. Evaluation of Appropriate Reference Genes for Reverse Transcription-Quantitative PCR Studies in Different Tissues of a Desert Poplar via Comparision of Different Algorithms. Int J Mol Sci. 2015;16(9):20468–91. pmid:26343648
- 49. Deguchi M, Bogush D, Weeden H, Spuhler Z, Potlakayala S, Kondo T, et al. Establishment and optimization of a hemp (Cannabis sativa L.) agroinfiltration system for gene expression and silencing studies. Scientific Reports. 2020;10(1):3504. pmid:32103049
- 50. Guerriero G, Mangeot-Peter L, Legay S, Behr M, Lutts S, Siddiqui KS, et al. Identification of fasciclin-like arabinogalactan proteins in textile hemp (Cannabis sativa L.): in silico analyses and gene expression patterns in different tissues. BMC Genomics. 2017;18(1):741-. pmid:28931375
- 51. Deguchi M, Kane S, Potlakayala S, George H, Proano R, Sheri V, et al. Metabolic Engineering Strategies of Industrial Hemp (Cannabis sativa L.): A Brief Review of the Advances and Challenges. 2020;11(1898).
- 52. Stout JM, Boubakir Z, Ambrose SJ, Purves RW, Page JE. The hexanoyl-CoA precursor for cannabinoid biosynthesis is formed by an acyl-activating enzyme in Cannabis sativa trichomes. 2012;71(3):353–65.