Tissue-specific evaluation of suitable reference genes for RT-qPCR in the pond snail, Lymnaea stagnalis

Care of snails

Animal use was consistent with the Canadian Council for Animal Care guidelines. A colony of L. stagnalis was bred and maintained in the animal care facility at StFX. The animals were exposed to a photoperiod matched to natural daylight patterns in Antigonish, Nova Scotia. The water in the animals’ tanks was changed three times per week. The animals were fed fish food and romaine lettuce ad libitum plus sinking protein pellets once per week.

Total RNA extraction

Total RNA was collected from each L. stagnalis tissue of interest. To prepare for RNA isolation, animals were anesthetized in 0.125% 1-phenoxy-2-propanol dissolved in Lymnaea saline for up to 30 min (Wyeth et al., 2009). The animals were dissected to isolate the six tissues of interest: CNS, tentacles, lips, penis, foot, and mantle. The tissues of each type from two animals were combined for each replicate to ensure adequate yields of RNA from each tissue, three replicates were produced with this method. Tissues were flash frozen in liquid nitrogen, shredded with razor blades and added to 500 µL of TRIzol reagent (Thermo Fisher, Waltham, MA, USA). Tissue solutions were thoroughly homogenized with a PowerGen 125 tissue homogenizer (Thermo Fisher, Waltham, MA, USA). Total RNA was extracted from the homogenized tissue via addition of 0.5 volumes of a 24:1 solution of chloroform and isoamyl alcohol. The aqueous layer was separated, and RNA was precipitated with isopropyl alcohol. RNA precipitate was transferred to the spin column of the E.Z.N.A.^® Mollusc RNA Kit (Omega Bio-Tek, Norcross, GA, USA) where the RNA was washed and subjected to the on-column DNAse I treatment as per the manufacturer’s instructions. To confirm the effectiveness of the DNase I treatment, PCR was performed with primer sets for Lst-tyrosine hydroxylase (forward: 5′-CCCACGTGTATCGTCATCTTA-3′; reverse: 5′-ATCTTCTCCTCGCAAAACC-3′) and Lst-EF1α with 20 ng of total RNA as a template as this reaction would not work in the absence of genomic DNA (gDNA) contamination. Lst-tyrosine hydroxylase was chosen as one of the gDNA contamination controls as it is frequently amplified by us in all of the tissues of interest as a part of a larger research project on neural development. No amplification was observed for either gene in any of the RNA samples after 40 cycles, thus the RNA samples were judged to be free of gDNA contamination.

Total RNA was quantified in a QuBit 3.0 (Thermo Fisher, Waltham, MA, USA). Yields of 50–200 ng/µL dissolved in 30 µL of molecular grade water were common for roughly 40 mg of tissue mass. A sample of total RNA was measured with a spectrophotometer to confirm a 260:280 ratio of 1.8–2.0. Lastly, five µL of each RNA sample was denatured in five µL of 2× RNA Loading Dye (New England Biolabs, Ipswich, MA, USA) at 75 °C for 10 min and then immediately placed on ice. The RNA samples were loaded on a 1% agarose non-denaturing gel and run alongside a single stranded RNA ladder (New England Biolabs, Ipswich, MA, USA) to verify that the 28S and 18S rRNA bands were intact and the RNA samples were not degraded.

Reverse transcription PCR

Reverse transcription-PCR was performed with iScript Reverse Transcription Supermix for RT-qPCR (Bio-Rad Laboratories, Hercules, CA, USA) which contains a mixture of oligo(dT) primers and random hexamers. For RT-PCR reactions, 500 ng of RNA was added to four µL of iScript Supermix and topped up to 20 µL with molecular grade water as recommended by the manufacturer. RT-PCR took place in a Bio-Rad CFX Connect thermocycler (Bio-Rad Laboratories, Hercules, CA, USA). The RT-PCR program ran at 25 °C for 5 min, 46 °C for 20 min and 95 °C for 1 min.

Design, validation, and optimization of primers for quantitative PCR

Candidate primers intended for qPCR were designed with CLC Main Workbench software (Qiagen, Hilden, Germany). Primers were designed to have a length of 19–23 nucleotides, a melting temperature between 55 and 62 °C, a GC content between 40% and 60% and produce a product no more than 200 bp long. Additionally, the GC content of each primer was kept within 10% of its counterpart, and the melting temperature of each was kept within 3 °C of its counterpart.

Two to six sets of primers were designed for each gene so that optimal primers could be selected for the qPCR reactions. The primers were tested on combined cDNA samples to minimize tissue-specific bias. Melt curves were performed to verify that one product was amplified. If the primers produced a single product, then five µL of the PCR reactions were run on a 2% (w/v) agarose gel at 60 V for roughly 30 min alongside a 50 bp DNA Ladder (New England BioLabs, Ipswich, MA, USA). The agarose gels were analyzed with a Bio-Rad ChemiDoc (Bio-Rad Laboratories, Hercules, CA, USA) to visualize the size of the product. Ultimately, one set of primers was selected for each candidate reference gene based on the quality of the PCR product. The primers for each candidate reference gene are listed in Table 1.

After the optimal primer sets were selected based on reaction specificity, primer efficiencies were calculated for each primer set. Primer efficiency curves were generated from RT-qPCR reactions on serial dilutions of RNA. The parent reaction contained 100 ng total RNA and four 1:5 dilutions were performed to generate five Cq values.

Quantitative PCR

The minimum information for publication of quantitative real-time (MIQE) PCR experiments guidelines were followed throughout the collection of qPCR data (Bustin et al., 2009). Amplification of all genes was detected with SyBR Green dye which generates fluorescence based on the synthesis of double-stranded DNA. The reactions contained two µL of cDNA with 10 µL of Bio-Rad SsoAdvanced Universal SyBR Mix, 600 nM forward and reverse primer concentration, and topped to 20 µL with DEPC H₂O. Each replicate of L. stagnalis tissue was subjected to qPCR reactions in triplicate. The qPCR reactions took place in a Bio-Rad CFX Connect thermocycler running a custom program. The custom qPCR program consisted of 95 °C for 30 s; 40 cycles of 95 °C for 15 s, 55 °C for 30 s. The plate was read by the machine to measure fluorescence at the end of each cycle.

Data analysis

The expression stability of each gene was assessed with four computational algorithms: geNorm (Vandesompele et al., 2002), NormFinder (Andersen, Jensen & Ørntoft, 2004), BestKeeper (Pfaffl et al., 2004), and RefFinder (Xie et al., 2012). GeNorm ranks candidates by their expression stability (M) values that are assigned following pair-wise variation measurements among genes. Genes with M < 1.5 are considered to be stable. The NormFinder algorithm compares intra-group variation (i.e., mRNA levels of one gene within the tissue of interest) to inter-group variation (i.e., mRNA levels of other genes in the same tissues) and assigns genes a stability value based on variation among Cq values. BestKeeper builds a stability index based on repeated pair-wise correlation analyses between every reference gene and judges reference gene stability based on the standard deviation (SD) from the index, genes that have a smaller SD will be ranked more highly by BestKeeper. Finally, RefFinder was used to combine the ranked results from each algorithm and assign each gene an overall rank.

Data analyses were performed to compare the stability of among candidate reference genes in each tissue of interest. Following the tissue-specific analyses, we also combined all data for analysis together to simulate the conditions of a whole-body analysis. This approach allowed for the identification of which genes are most stable in each tissue of interest and provided evidence for which genes would be the strongest reference genes in an RT-qPCR experiment that used RNA from a whole L. stagnalis body. Given that a minimum of two reference genes are typically recommended for proper RT-qPCR data normalization, the recommended reference genes are generally expressed in terms of pairs.

Primer specificity and efficiency

Primers for all candidate reference genes were evaluated to ensure that they could produce consistent results and not amplify off-target products or generate primer dimers. Following amplification, each primer pair produced amplicons that yielded single bands at the correct size after electrophoresis in 2% agarose gels (Fig. 1). Additionally, no amplification was observed in controls that lacked reverse transcriptase in the RT-PCR or lacked cDNA template or primers in qPCR. Thus, primer pairs specifically amplified a single cDNA target. Based on the standard curves, primer set efficiencies ranged from 90% (GAPDH and Rapac1) to 106% (ACTB) with correlation coefficients (R²) of >0.980 (Table 1).

Tissue-specific expression profiles and stability of reference genes

The candidate reference genes showed a variety of expression levels across individual tissues (Fig. 2). Lst-UBI produced the smallest mean Cq values and displayed the smallest Cq ranges for most tissues (Table 2). Conversely, Lst-VGKC mRNA was the least abundant in every tissue and displayed large Cq ranges.

The algorithms used to measure the stability of each candidate reference gene were highly congruent in their tissue-specific rankings. Therefore, only the results from the RefFinder analysis (which combines the rankings from the other three algorithms) are presented in Fig. 3 (results for other algorithms reported in Figs. S1–S3). Lst-GAPDH and Lst-EF1α were the most stable candidate reference genes in the CNS and penis. Lst-UBI and Lst-GAPDH were the most stable in the tentacles, Lst-GAPDH and Lst-Rapac1 were the most stable in the lips, Lst-UBI and Lst-ACTB were the most stable in the foot, and Lst-UBI and Lst-VGKC were the most stable in the mantle. Overall, identification of the most reliable candidate reference gene was tissue dependent. However, according to geNorm, all genes except for Lst-ACTB and Lst-TUBB were acceptable as a reference genes for all tissues as they met the minimum stability threshold of M < 1.50 for heterogeneous tissue samples.

Expression levels and stability of reference genes in a combined tissue analysis

To assess which candidate reference genes could be the most stable in whole-body preparations of L. stagnalis, the data from all tissues was combined and analyzed together with geNorm, NormFinder, BestKeeper, and RefFinder. The candidate reference genes spanned a wide range of expression levels when data from all tissue types was combined. Based on the Cq values, Lst-UBI was the most highly expressed gene (mean Cq 16.62 ± 0.09) and had the smallest overall Cq range whereas Lst-VGKC was the least expressed gene and had the highest biological variability as measured by the standard error of the mean (mean Cq 29.87 ± 0.19; Table 2). All four algorithms used to assess expression stability were highly congruent in their ranking of candidate reference genes in the combined analysis (Fig. 4). GeNorm, NormFinder, and RefFinder rated Lst-GAPDH and Lst-EF1α as the top two most stable reference genes. BestKeeper ranked Lst-UBI as the most stable with Lst-EF1α and Lst-GAPDH as the second and third most stable, respectively.