Accepted Manuscript Validation of housekeeping genes for quantitative mRNA expression analysis in OsHV-1 infected ark clam, Scapharca broughtonii Lusheng Xin, Bowen Huang, Changming Bai, Chongming Wang PII: DOI: Reference:
S0022-2011(17)30436-6 https://doi.org/10.1016/j.jip.2018.04.011 YJIPA 7093
To appear in:
Journal of Invertebrate Pathology
Received Date: Revised Date: Accepted Date:
19 October 2017 24 April 2018 30 April 2018
Please cite this article as: Xin, L., Huang, B., Bai, C., Wang, C., Validation of housekeeping genes for quantitative mRNA expression analysis in OsHV-1 infected ark clam, Scapharca broughtonii, Journal of Invertebrate Pathology (2018), doi: https://doi.org/10.1016/j.jip.2018.04.011
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Validation of housekeeping genes for quantitative mRNA expression analysis in OsHV-1 infected ark clam, Scapharca broughtonii Lusheng Xin a, b,1, Bowen Huang c,1, Changming Bai a, b, Chongming Wang a, b* a
Qingdao Key Laboratory of Mariculture Epidemiology and Biosecurity, Key
Laboratory of Maricultural Organism Disease Control, Ministry of Agriculture, Yellow Sea Fisheries Research Institute, Chinese Academy of Fishery Sciences, Qingdao 266071, China; b
Function Laboratory for Marine Fisheries Science and Food Production Processes,
Qingdao National Laboratory for Marine Science and Technology, Qingdao 266071, China c
Tianjing Agriculture University, Tianjing 300191, China
*
Correspondence authors:
Tel: 86-532-85823062 E-mail:
[email protected] (C. Wang) 1
These authors contributed equally to this work.
1
Abstract Ostreid herpesvirus-1 (OsHV-1) presents interspecies transmission among bivalves. Recently, events of mass mortalities of ark clams (Scapharca broughtonii) infected with OsHV-1 have been recorded. To accurately assess the gene responding patterns of ark clams post OsHV-1 infection, constant stable housekeeping genes (HKGs) are needed as internal control to normalize raw mRNA expression data. In this study, ten candidate HKGs were selected, including 18S rRNA (18S), beta-actin (ACT), Glyceraldehyde 3-phosphate dehydrogenase (GAPDH), NADH dehydrogenase subunit (NADH), Elongation factor-1a (EF-1a), Elongation factor-1β (EF-1β), Elongation factor-1γ (EF-1γ), Ribosomal protein L7 (RL7), Ribosomal protein L15 (RL15) and Ribosomal protein S18 (S18). The expression levels of ten candidate HKGs were analyzed by real-time PCR under given experimental conditions, including various tissues, OsHV-1 challenge, temperature stress and OsHV-1 challenge at different temperature. Their expression stability values were further calculated using two different statistical models (geNorm and NormFinder). The results showed that different tissues presented distinct best pair genes combinations for gene expression analysis under OsHV-1 challenge. RL15 was comparatively more stable than other HKGs under various experimental conditions, while commonly used 18s and ACT seemed to be more greatly influenced by most given experimental conditions in ark clams. This study emphasized the necessity of prior validation of HKGs and would facilitate future gene expression analysis in ark clams or other shellfishes. Keywords: housekeeping gene; Ostreid Herpesvirus; Scapharca broughtonii; shellfish; 2
Real-time PCR
1. Introduction Ostreid herpesvirus 1 (OsHV-1) was characterized as a major class of the herpesviruses (Davison et al., 2005), and emerged as a disastrous pathogen in shellfish aquaculture. OsHV-1 was often involved with mass mortality events of ostreid species (Keeling et al., 2014; Segarra et al., 2010). Besides vertical transmission across offspring (Barbosasolomieu et al., 2005), OsHV-1 appears capable of interspecies transmission (Arzul et al., 2001). Subsequently, other shellfish species were also characterized with OsHV-1 infection cases, including ark clam (Scapharca broughtonii) (Xia et al., 2015). Ark clam is distributed widely along the northwest pacific coast, and has become one of the most commercially important shellfish cultivated in China and other northern Asia countries (Yokogawa, 1997). Mass mortalities of ark clams were recorded and resulted in great economic loss during the early summer of 2012 and 2013, and had been confirmed to be associated with OsHV-1 infection (Bai et al., 2016). Until now, the interaction between pathogen OsHV-1 and host ark clam is largely unknown, specific responding genes as biomarkers remain to be further explored for early warning and effective management against OsHV-1 prevalence in ark clams. Though several methods have been established for the diagnosis of OsHV-1 infection (e.g., transmission electron microscopy, in situ hybridization and PCR diagnostic methods) (Barbosa-Solomieu et al., 2004; Oden et al., 2011; Pepin et al., 2008), these methods could only sense OsHV-1 particles after onset and are not available for 3
distinguishing vulnerable population. Specifically responding genes in OsHV-1 infected ark clams are expected to become much earlier, more sensitive and effective diagnostic objects for OsHV-1 detection, and further guide the acknowledgement of host-pathogen interaction. For searching these genes, mRNA expression profile before and after OsHV-1 infection is a critical cue, which could be assessed by semi-quantitative PCR, Northern Blot, RT-PCR and so on. In comparison, RT-PCR is more frequently used for its high convenience, sensitivity and reproducibility (Nolan et al., 2006). Even mRNA expression levels were quantified by costly high-throughput sequencing, such as microarray and RNA-seq, still further need to be confirmed by RT-PCR. Thus RT-PCR serves an important tool for identifying genes responding to OsHV-1 infection. Due to the unmanageable variations of initial total RNA amount and cDNA synthesis efficiency, constant stable expression HKGs are critical to normalize raw data for RT-PCR analysis (Dheda et al., 2004). Microsoft Excel geNorm and NormFinder were mostly employed to determine HKG stability and suitability under given experimental conditions, ranking the candidate genes according to the calculated expression stability values, more stable genes present lower expression stability values (Gutierrez et al., 2008; Jensen and Ørntoft, 2004; Pan et al., 2015). The accuracy and repeatability of target gene expression results are strongly dependent on selected HKGs. Least stability HKGs could even cover the significant fluctuation of target genes expression (Lacroix et al., 2014). However, limited previous studies were reported about the selection of HKGs in molluscs (Table. 2). 4
Available HKGs were only in several shellfish species, such as Crassostrea gigas (Dheilly et al., 2011), Ostrea edulis (Morga et al., 2010) and Mytilus edulis (Cubero-Leon et al., 2011). In addition, to date, most validated HKGs in molluscs exposed to pathogens (Vibrio splendidus, Bonamia ostreae) were tested in haemocytes only (Araya et al., 2008; Morga et al., 2010). HKGs were validated in oyster larvae and spat under OsHV-1 infection (Du et al., 2013; Segarra et al., 2014). Many of these HKGs may not be stably expressed across species under different conditions (Araya et al., 2008; Mauriz et al., 2012). Thus prior validation is required when a HKG is suggested to be used in a new given experimental condition and can not be ignored. There are no validated HKGs that have been suggested in ark clams up to now. Whereas previous studies mainly selected HKGs, β-actin or 28s, for gene expression analyses in ark clams without validation (Zheng et al., 2016). The suitability of these two HKGs was doubtful, because they exhibited significant fluctuation of mRNA expression levels and were abandoned under some given experimental conditions, such as white spot syndrome virus (WSSV) challenge shrimp and mussel various stages of gametogenesis (Cubero-Leon et al., 2011; Dhar et al., 2009). In this study, 10 candidate HKGs were selected based on molecular function and stable gene expression levels from our previous RNA-seq data. Their expression stability values were further assessed among different tissues and OsHV-1 challenge. Temperature, an important trigger factor of OsHV-1 disease outbreak (above 13°C) (Pernet et al., 2015; Petton et al., 2013), was also taken into account as the given experimental condition under that HKG expression stability was checked. The aim of this study was to 5
suggest available HKGs, referring to tissue distribution, OsHV-1 challenge or temperature stress in ark clams.
2. Materials and Method 2.1 Ark clams and tissue sampling Healthy ark clams (two-year-old) with an average shell length of 55mm were obtained from a local farm in Qingdao, China. Prior to the following treatment, ark clams were acclimated in tanks (10 ark shells in each tank) with aerated sea water at 18.0°C for 13 days. For tissue distribution analysis of gene expression, tissue samples (including haemocyte, mantle, gill, digestive gland, adductor muscle and foot) were collected from 6 ark clams and placed into TRIzol (Invitrogen, USA) for later RNA isolation. Three replicates were employed for each tissue sample. 2.2 Experimental infection and temperature stress For experimental infection, viral inoculum was firstly prepared using OsHV-1 infected ark clams as previous description using oyster with modification (Schikorski et al., 2011). Briefly, mantle tissue (50mg each individual) was sampled from six OsHV-1 infected ark clams and homogenized. After a short centrifugation (4°C, 1000g), the supernatant of mantle homogenate was remained and successively filtered through a series of syringe filters with pore size at 5μm, 2μm, 0.45μm and 0.22μm. The filtered tissue homogenate was used as a viral inoculum and stored for 2 days at 4°C. Meanwhile, un-infected ark clams were used to make negative tissue 6
homogenates as the negative control. After that, 200 μL of each tissue homogenate was used for DNA extraction and OsHV-1 quantification by qPCR according to the previously described protocol (Bai et al., 2016). Two hundred and forty ark clams were randomly selected and separated into two groups (HN and HP). Both were remained at 18°C, group HP was injected with 100 μL viral inoculum (~105 copies of viral DNA/μL) into the foot of each ark clam. Correspondingly the other group HN was injected with the same volume of negative tissue homogenate as the negative control. Tissue samples (haemocyte, mantle, gill, digestive gland, adductor muscle and foot) of six ark clams were respectively collected from each group at 0h and every 12h for 72h post injection. For temperature stress analysis, after 13 days at 18 °C, one hundred arks shells were transferred to aerated sea water at 12°C. After continuous culture for 48h, haemocyte and mantle samples were collected from six ark clams at 12°C and parallel six at 18°C respectively. In order to analyze gene expression under OsHV-1 challenge at different temperature, ninety ark clams were also randomly selected and separated into two groups (LN and LP) before being transferred to aerated sea water at 12°C, group LP received the same injection as group HP, while group LN as group HN. After undergoing 48h culture at 12°C, haemocyte and mantle samples of six ark clams were also collected from group LN and LP, to compare with group HN and HP remained at 18°C for 48h . All collected tissue samples were paralleled with three replicates. 2.3 RNA extraction and cDNA synthesis 7
Total RNA was extracted from tissue samples using TRIzol reagent (Invitrogen, USA) according to the manufacturer’s protocol. Genomic DNA was eliminated by digestion with RNase-free DNase I (Promega). The quality and quantity of RNA samples were assessed by Nanodrop 2000 spectrophotometer (Thermo). Only RNA samples with a 260/280 ratio among 1.9-2.1 and a 260/230 ratio among 2.0-2.5 were used for analysis. The integrity of RNA samples was confirmed by agarose gel electrophoresis. First-strand cDNA was synthesized based on M-MLV reverse transcriptase (Promega) by using 1 μg total RNA as template according to the manufacturer’s instructions. The cDNA products were stored at -80°C for the subsequent real-time PCR. 2.4 Selection of candidate HKGs and primer design Ten candidate HKGs were selected according to their functions and expression fluctuations in our previous transcriptome data of ark clam haemocytes under given experimental conditions. Nucleotide sequences of ACT, GAPDH, NADH, EF-1a, EF-1β, EF-1γ, RL7, RL15 and S18 were obtained from our previous transcriptome data and confirmed by cloning and sequencing, while that of 18s was gained from GenBank database (GenBank:JN974500.1). Primers were designed using Primer Premier 5 software (PREMIER Biosoft). All the amplicons of selected candidate HKGs had been purified from agarose gel, and ligated to pMD 19-T Vector (TAKARA) for sequencing with six replicates in Sangon Biotech Company, China. Detailed information on these candidate HKGs and all the primers used for real-time PCR were shown in Table 1. 8
2.5 Real-time PCR analysis Real-time PCR was performed using Bio-Rad CFX Connect RealTime system (Bio-Rad Laboratories). Firstly, primer specificities were checked using melting curves and gel picture of corresponding amplicons. PCR efficiency was calculated for each primer pair by making standard curves of 10-fold diluted cDNA (1 to 1:10,000), using the formula: E = 10-1/slope (Pfaffl, 2001). A 20μl RT-PCR reaction volume contained 2 μL diluted cDNA template (1:10), 10 μL 2×SYBR Green PCR Mix (Takara), 0.5 μL of each gene-specific primer (250 nM) and 7 μL distilled water. Cycling conditions were as follows: 95°C for 30s, 40 cycles of 95°C for 5 s and 65°C for 35s and then a melt curve stage: 20 s at 95 °C and 30 s at 60 °C followed by temperature increasing (0.5°C per 10s) into 95°C. RT-PCR was performed in triplicate for each sample. 2.6 Gene expression stability analysis The Cq values were transformed to linear scale expression quantities as the form of 2−ΔCq. ΔCq value equals Cq value of HKGs obtained from each sample subtracting the lowest one (ΔCq= Cqeach- Cqmin). Then the data that stands for HKG expression level was input into and analyzed by Microsoft Excel geNorm (version win_3.5) (Pattyn et al., 2003) and NormFinder (version 0953) (Jensen and Ørntoft, 2004). Gene expression stability values were directly output. Genes with the lowest stability value have the most stable expression. Statistical analyses were conducted with the IBM SPSS Statistics 19.0 package, the box-and-whisker plots were obtained by this software. 9
3. Result 3.1. Selected candidate HKGs Ten candidate HKGs, including ACT, GAPDH, NADH, EF-1a, EF-1β, EF-1γ, RL7, RL15, S18 and 18s, were selected and further assessed. NADH, EF-1a, EF-1β, EF-1γ, RL7, RL15, and S18 showed parallel colors in the heatmap (Fig.1), indicating similar expression profiles under given experimental conditions based on our previous RNA-seq data analysis. Though ACT and GAPDH emerged distinct color difference in the heatmap (Fig.1), their molecular functions were fundamental and indispensable to sustain cells, and had been widely used as HKGs together with 18s. Thus ACT, GAPDH and 18s were also taken into account in our studies. Besides nucleotide sequence of 18s had been released and directly gained in the GeneBank database (accession: JN974499), the others were acquired from our RNA-seq assembled unigenes. All selected candidate HKGs exhibited high identities (>70% in nucleotide sequence) to their respective homologs in NCBI database, the sequences of ten candidate HKGs were cloned and confirmed by sequencing to facilitate subsequent RT-PCR primer design. 3.2. Amplification specificity and average expression levels of candidate HKGs The specificity of RT-PCR primers was confirmed by a signal band in each lane of the agarose gel with expected amplicon length (117-383bp) (Fig.2), and by corresponding melting curve with a single peak (Fig.S2). All the amplicons of candidate HKGs were confirmed by sequencing. Amplification efficiencies of 10
RT-PCR primers ranged between 95.9% and 104.6% (Table 3), S18 showed the lowest (95.9 %). The Pearson's moment correlation coefficients were 0.990. Box plots showed the comparative expression levels of ten candidate HKGs (Cq values) among six tissues (haemocyte, mantle, gill, digestive gland, adductor muscles and foot) in ark clams (Fig.3). ACT and 18s exhibited higher degree of dispersion with greater box lengths among ten candidate HKGs. Meanwhile, 18s showed the lowest Cq median indicating the highest average gene expression level among six tissues, which was in accordance with the lowest Cq value of 18s measured in the mixed six tissue samples (Fig.S1). 3.3. Expression stability of candidate HKGs among various tissues Among various tissues, geNorm was further used to calculate the gene expression stability values (M) as shown in Fig.4. The maximum threshold of M value was recommended at 1.5 for suitable HKGs selection, lower M value presents greater gene stability. Among various tissues, the M values of 18s, ACT, EF-1β and S18 (M=3.07, 2.49, 1.74 and 1.55, respectively) were above the threshold (shown by the dot line), while the other six were below that (Fig.4a). Gene pairs RL7 and EF-1a emerged as the most stable genes across tissues. Gene expression stability was also evaluated by NormFinder based on their minimum variance combined intra and inter groups (model-based approach) as shown in Fig.7, more stable genes also present lower values of average expression stability (listed in Table.3). NormFinder also ranked 18s, ACT, EF-1β and S18 as the top four fluctuant genes among different tissues, and RL7 represented in the top three stable genes. 11
To determine the optimal number of reference genes required for normalization, geNorm software was used by means of pairwise variation (Vn/n+1) calculated between the normalization factors NFn and NFn+1( using n or n+1 reference genes respectively)(Vandesompele et al., 2002). A cut-off value for the pairwise variation was suggested at 0.15, below this value indicating an additional HKG was not required. For example, V2/3 <0.15, two HKGs would satisfy the analysis; V2/3 >0.15 and V3/4<0.15, three HKGs were needed. However, among various tissues of ark clams, the values of pairwise variation Vn/n+1 were all out of range (the least two, V3/4 =0.230 and V5/6=0.230) (Fig.4b). According to the geNorm software manual, the proposed 0.15 value must not be taken as a too strict cut-off. Three best HKGs (RL7, EF-1a and EF-1γ) were suggested for gene expression analysis in this case. 3.4. Expression stability of candidate HKGs under stress Under the given condition of OsHV-1 challenge, the M values of ten candidate HKGs were calculated by geNorm (Fig.5). In haemocytes, mantle tissue, adductor muscles and foots, pair genes RL7 and RL15 (M=0.206), S18 and EF-1β (M=0.192), EF-1a and EF-1γ (M=0.224), S18 and EF-1β (M=0.154) respectively emerged as the most stable combinations under OsHV-1 challenge (Fig.5A, B, C and D). The values of pairwise variation V2/3 were all below the cut-off value (0.15) in above four tissues under OsHV-1 challenge (Fig.5a, b, c and d), which indicated two genes could satisfy RT-PCR data normalization. Under temperature stress, the M values of ten candidate HKGs were ranged from 0.044 to 0.376 and 0.061 to 0.233 in haemocytes and mantle tissue, respectively (Fig.6 A and B), and ranged from 0.060 to 0.413 and 12
0.040 to 0.242 respectively under OsHV-1 challenge at different temperature (Fig.6 C and D). The values of pairwise variation V2/3 were also all below 0.15 in haemocytes and mantle tissue under above two conditions (Fig.6 a, b, c and d), indicating two HKGs could satisfy gene expression analyses. In gills and digestive glands, the mRNA expression levels of three candidate HKGs (ACT, NADH and GAPDH) were not detected at some time points under OsHV-1 challenge and these genes were not taken into account for the subsequent analyses in gills and digestive glands (Fig.5 E and F).The M values of the remained seven candidate HKGs were ranked in gills as follows: 18s (M=1.651), RL7 (M=1.250), EF-1γ (M=1.200), RL15 (M=1.120), EF-1a (M=0.988), S18 (M=0.473) and EF-1β (M=0.473) (Fig.5 E); while in digestive gland as follows: RL7 (M=1.591), RL15 (M=1.469), 18s (M=1.131), EF-1a (M=0.871), EF-1β (M=0.685), EF-1γ (M=0.318) and S18 (M=0.318) (Fig.5 F). The values of pairwise variation Vn/n+1 were all above 0.15 in gills and digestive glands. V5/6 (0.186) and V3/4 (0.243) appeared the lowest in gills and digestive glands respectively (Fig.5 e and f). The expression stability values of candidate HKGs under stress were also calculated by NormFinder (Table.3), and the ranking of ten candidate HKGs under given experimental conditions were similar to the results calculated by geNorm. S18 and EF-1γ presented most frequently (six and seven times, respectively) in top three stable genes in tissue samples under stress, flowed by EF-1a, EF-1β and RL15 that all presented in top three for four times (shown in Table.3). Under tested experimental conditions, 18s and ACT showed much larger amount of expression stability values 13
(total up to 7.215 and 5.034) than the others, while RL15 was the smallest (total up to 1.646) (Fig.7).
4. Discussion The current study selected ten candidate HKGs from ark clams, and assessed the influence of factors, tissue distribution, OsHV-1 challenge and temperature stress, on the expression stability of ten candidate HKGs. In ark clams, pair genes RL7 and EF-1a were assessed as the most stable genes among various tissues. While pair genes GAPDH and ARF1 (ADP-ribosylation factor 1) in oyster, pair genes RPL5 (ribosomal protein L5) and CY (cyclophilin) in abalone appeared most stable among various tisssues (Table. 2). Besides not being assessed in oysters, 60S ribosomal protein subunit related genes (RL7 and RPL5) seemed more stable across tissues than others in mollusk species. Some common HKGs gene expression levels exhibited great variation among different tissues, especially 18s and ACT which were least stable (based on stability values) in ark clams. Though 18s and ACT are commonly used as internal controls for gene expression analysis (Deane and Woo, 2005; Kang et al., 2006), their suitability and stability should be re-assessed. To our knowledge, our research was the first to determine the suitability of the 18S gene across tissue types in mollusk species. ACT presented as one of top two least stable genes among various tissues in oyster and abalone (Dheilly et al., 2011; López-Landavery et al., 2014), which was different from fish (such as zeabrafish and Japanese flounder) ACT exhibiting more stable across tissues among tested candidate 14
HKGs (Mccurley and Callard, 2008; Zheng and Sun, 2011). HKGs referring to pathogen challenging have been selected in several mollusc species, whereas mainly in haemocytes (Table. 2), related RT-PCR primer information of validated HKGs was listed in Table.S1. EF-1a and GAPDH were characterized as the most stable genes in haemocytes of Ostrea edulis under parasite Bonamia ostreae challenging in vitro (Morga et al., 2010), while EF-1a and S18 exhibited the most stability in haemocytes of Mya arenaria under bacteria Vibrio splendidus challenging in vitro (Araya et al., 2008). Selection of HKGs referring to OsHV-1 challenge was also reported in oyster larva and spat mantle, which suggested RL7 and S18 as the most suitable ones in oyster larva (Du et al., 2013), EF-1a in oyster spat mantle(Segarra et al., 2014). Thus these HKGs were all compared in our study on ark clams post OsHV-1 challenge, different tissues presented distinct best pair genes combination (Fig.5). In haemocytes, mantle tissue, foots and adductor muscles of ark clams, the best pair genes could satisfy gene expression analysis post OsHV-1 challenge, S18 and EF-1β presented more stable than EF-1a in ark clam mantle that were not tested in oyster spat mantle(Segarra et al., 2014). However, in gills and digestive glands of ark clams, gene expression levels might be greatly influenced by OsHV-1 challenge, most of tested candidate HKGs appeared significant fluctuations. According to the results of pairwise variation Vn/n+1 calculated by geNorm, more than two HKGs were recommended to be used for RT-PCR analysis in gills and digestive glands post OsHV-1 challenge. What’s more, all the most stable genes in different tissues of ark clams post OsHV-1 challenge belong to two families, 15
ribosomal protein subunit genes and elongation factor genes. It was almost the same as validated HKGs referring to pathogen challenging in other shellfish species (Araya et al., 2008; Du et al., 2013; Morga et al., 2010). In another important aquaculture crustacean shrimp (Penaeus stylirostris) under virus (white spot syndrome virus, WSSV) challenging, elongation factor gene EF-1a also appeared as the most stable HKG (Dhar et al., 2009). Therefore, ribosomal protein subunit genes and elongation factor genes take priority to be recommend for the selection of HKGs under virus or other pathogen challenging in molluscs. Different temperatures (12°C and 18°C) showed less influence on the expression of ten candidate HKGs in ark clams’ haemocytes and mantle tissue. The stability values of ten candidate HKGs were all far below the threshold (M=1.5), which indicated that ten candidate HKGs were all eligible for gene expression analysis under this condition. Under OsHV-1 challenge at different temperature (12°C and 18°C), pair genes RL15 and NADH, EF-1γ and S18, were recommended as the most stable genes in ark clams’ haemocytes and mantle tissue, respectively, which is similar to under OsHV-1 challenge only. OsHV-1 challenge took up the dominant factor that affected HKGs’ expression when ark clams were under OsHV-1 challenge at different temperature (12°C and 18°C). In conclusion, OsHV-1 has gained comprehensive attention as a highly pathogenic virus which transmits between and among marine shellfish species causing mass mortalities (Arzul et al., 2001; Keeling et al., 2014; Schikorski et al., 2011; Xia et al., 2015). Characterization of specific responding genes of ark clams post OsHV-1 16
infection or temperature stress would provide a new approach for early warning and detecting of OsHV-1 infection, and enrich our knowledge about the interaction between OsHV-1 pathogens and shellfish hosts. Prior selection of HKGs were critical for gene expression analysis to identify these specific responding genes. Ten candidate HKGs of ark clams were investigated, and suitable HKGs were recommended referring to different tissues, OsHV-1 challenge or temperature stress. Considering all the test conditions, commonly used HKGs 18s and ACT gained the two highest total stability values, seemed to be more greatly influenced by most given experimental conditions. 18s and ACT as HKGs were also criticized in other molluscs (Morga et al., 2010; Siah et al., 2008), should be cautiously selected for gene expression analysis. RL15 exhibited the lowest, which suggested RL15 was comparatively more stable than other HKGs under various experimental conditions, could be firstly selected and assessed in ark clams or other molluscs under a new given experimental conditions. Except the factor of different tissues, S18 and EF-1γ appeared most frequently in top three most stable genes under tested given conditions, exhibited more stable expression levels and were also recommended for further studies of HKGs selection in an individual tissue of molluscs. The results of this study would facilitate gene expression analysis in ark clams and were helpful for further studies of HKGs selection in molluscs.
Acknowledgements We are grateful to all the laboratory members for their technical advice and 17
helpful discussions. This research was supported by the earmarked fund (CARS-48) for Modern Agro-industry Technology Research System, fund from the Special Scientific Research Funds for Central Non-profit Institutes, Yellow Sea Fisheries Research Institutes (Project No. 20603022017007).
Reference Araya, M. T., et al., 2008. Selection and evaluation of housekeeping genes for haemocytes of soft-shell clams (Mya arenaria) challenged with Vibrio splendidus. Journal of invertebrate pathology. 99, 326-331. Arzul, I., et al., 2001. French scallops: a new host for ostreid herpesvirus-1. Virology. 290, 342-9. Bai, C., et al., 2016. Identification and characterization of ostreid herpesvirus 1 associated with massive mortalities of Scapharca broughtonii broodstocks in China. Diseases of Aquatic Organisms. 118, 65-75. Barbosa-Solomieu, V., et al., 2004. Diagnosis of Ostreid herpesvirus 1 in fixed paraffin-embedded archival samples using PCR and in situ hybridisation. Journal of virological methods. 119, 65-72. Barbosasolomieu, V., et al., 2005. Ostreid herpesvirus 1 (OsHV-1) detection among three successive generations of Pacific oysters (Crassostrea gigas). Virus Research. 107, 47-56. Cubero-Leon, E., et al., 2011. Reference gene selection for qPCR in mussel, Mytilus edulis, during gametogenesis and exogenous estrogen exposure. Environmental Science & Pollution Research. 19, 2728. Davison, A. J., et al., 2005. A novel class of herpesvirus with bivalve hosts. Journal of General Virology. 86, 41-53. Deane, E. E., Woo, N. Y. S., 2005. Cloning and characterization of the hsp70 multigene family from silver sea bream: Modulated gene expression between warm and cold temperature acclimation. Biochemical and Biophysical Research Communications. 330, 776-783. Dhar, A. K., et al., 2009. Validation of reference genes for quantitative measurement of immune gene expression in shrimp. Molecular Immunology. 46, 1688-1695. Dheda, K., et al., 2004. Validation of housekeeping genes for normalizing RNA expression in real-time PCR. Biotechniques. 37, 118-9. Dheilly, N. M., et al., 2011. Development of a Pacific oyster (Crassostrea gigas) 31,918-feature microarray: identification of reference genes and tissue-enriched expression patterns. BMC genomics. 12, 468. Du, Y., et al., 2013. Validation of housekeeping genes as internal controls for studying gene expression during Pacific oyster (Crassostrea gigas) development by quantitative real-time PCR. Fish & Shellfish Immunology. 34, 939. Gutierrez, L., et al., 2008. Towards a Systematic Validation of References in Real-Time RT-PCR. The Plant Cell. 20, 1734-1735. Jensen, J., Ørntoft, T., 2004. Normalization of real-time quantitative RT-PCR data: a model based 18
variance estimation approach to identify genes suited for normalization-applied to bladder-and colon-cancer data-sets. Cancer Research. 64, 5245-5250. Kang, Y.-S., et al., 2006. Analysis of EST and lectin expressions in hemocytes of Manila clams (Ruditapes philippinarum) (Bivalvia: Mollusca) infected with Perkinsus olseni. Developmental & Comparative Immunology. 30, 1119-1131. Keeling, S. E., et al., 2014. New Zealand juvenile oyster mortality associated with ostreid herpesvirus 1 an opportunistic longitudinal study. Diseases of Aquatic Organisms. 109, 231-239. López-Landavery, E. A., et al., 2014. Selection of reference genes as internal controls for gene expression in tissues of red abalone Haliotis rufescens (Mollusca, Vetigastropoda; Swainson, 1822). Gene. 549, 258-65. Lacroix, C., et al., 2014. A selection of reference genes and early-warning mRNA biomarkers for environmental monitoring using Mytilus spp. as sentinel species. Marine pollution bulletin. 86, 304-313. Mauriz, O., et al., 2012. Selection of reference genes for quantitative RT-PCR studies on the gonad of the bivalve mollusc Pecten maximus L. Aquaculture. 370, 158-165. Mccurley, A. T., Callard, G. V., 2008. Characterization of housekeeping genes in zebrafish: male-female differences and effects of tissue type, developmental stage and chemical treatment. BMC Molecular Biology,9,1(2008-11-12). 9, 102. Morga, B., et al., 2010. Identification of genes from flat oyster Ostrea edulis as suitable housekeeping genes for quantitative real time PCR. Fish & shellfish immunology. 29, 937-945. Nolan, T., et al., 2006. Quantification of mRNA using real-time RT-PCR. Nature protocols. 1, 1559. Oden, E., et al., 2011. Quantification of ostreid herpesvirus 1 (OsHV-1) in Crassostrea gigas by real-time PCR: determination of a viral load threshold to prevent summer mortalities. Aquaculture. 317, 27-31. Pan, H., et al., 2015. A Comprehensive Selection of Reference Genes for RT-qPCR Analysis in a Predatory Lady Beetle, Hippodamia convergens (Coleoptera: Coccinellidae). PLOS ONE. 10, e0125868. Pattyn, F., et al., 2003. RTPrimerDB: the Real-Time PCR primer and probe database. Nucleic Acids Research. 31, 122-123. Pepin, J.-F., et al., 2008. Rapid and sensitive detection of ostreid herpesvirus 1 in oyster samples by real-time PCR. Journal of Virological Methods. 149, 269-276. Pernet, F., et al., 2015. Influence of low temperatures on the survival of the Pacific oyster ( Crassostrea gigas ) infected with Ostreid herpesvirus type 1. Aquaculture. 445, 57-62. Petton, B., et al., 2013. Petton B, Pernet F, Robert R, Boudry P.. Temperature influence on pathogen transmission and subsequent mortalities in juvenile Pacific oysters Crassostrea gigas. Aquacult Environ Interact 3: 257-273. Aquaculture Environment Interactions. 3, 257-273. Pfaffl, M. W., 2001. A new mathematical model for relative quantification in real-time RT–PCR. Nucleic Acids Research. 29, e45-e45. Schikorski, D., et al., 2011. Experimental infection of Pacific oyster Crassostrea gigas spat by ostreid herpesvirus 1: demonstration of oyster spat susceptibility. Veterinary Research. 42, 27. Segarra, A., et al., 2014. Transcriptomic study of 39 ostreid herpesvirus 1 genes during an experimental infection. Journal of Invertebrate Pathology. 119, 5-11. Segarra, A., et al., 2010. Detection and description of a particular Ostreid herpesvirus 1 genotype associated with massive mortality outbreaks of Pacific oysters, Crassostrea gigas, in France in 19
2008. Virus Research. 153, 92-99. Siah, A., et al., 2008. Selecting a set of housekeeping genes for quantitative real-time PCR in normal and tetraploid haemocytes of soft-shell clams, Mya arenaria. Fish & Shellfish Immunology. 25, 202-207. Vandesompele, J., et al., 2002. Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. Genome Biology. 3. Xia, J., et al., 2015. Complete genome sequence of Ostreid herpesvirus-1 associated with mortalities of Scapharca broughtonii broodstocks. Virology Journal. 12, 110. Yokogawa, K., 1997. Morphological and genetic differences between Japanese and Chinese red ark shell Scapharca broughtonii. Fisheries science. 63, 332-337. Zheng, L., et al., 2016. Ferritin has an important immune function in the ark shell Scapharca broughtonii. Developmental & Comparative Immunology. 59, 15-24. Zheng, W. J., Sun, L., 2011. Evaluation of housekeeping genes as references for quantitative real time RT-PCR analysis of gene expression in Japanese flounder (Paralichthys olivaceus). Fish & Shellfish Immunology. 30, 638.
20
Table 1 Quantity real time PCR primers’ information for candidate HKGs Gene symbol Actin 18s GAPDH NADH
EF-1α
EF-1β
EF-1γ
RL7
RL15
S18
Molecular function
Cytoskeletal structural protein
Primer equences (5’-3’)
F: GGTTGTGAATGAATAACCACGTTC R: ACCTTCTACAACGAACTCCGTGTC
Ribosomal RNA
F: CTCTTTGGTGGTTACTGCTCGTC R: TCGCTAGTTGGCATCGTTTATG
Glycolytic enzyme
F: AACAGCGACCACTTCTACATC R: TGAAATCCTTTATCCAGCAAC
NADH dehydrogenase subunit
F: ATTACCCTGAAGGTCACAAG R: CACTCCTCCAAGCATACAAG
translation elongation factor
F: TTAGACGCAATCCTCCCTCC R: TAACGACGGTTCCTGGTTTG
translation elongation factor
F: AATTGTCCTTGATGTGAAGCCTTG R: GTTAAATGCAACAACGTCGACTGA
translation elongation factor
F: GTATGTGCGTGTTTGTAGGTACGT R: ATTATGTTGGCAATGCACAGTTAC
60S ribosomal protein subunit
F: AGATTTGTAAGATGGCCCAAGTAC R: CAAGAGGTTCAACATCATGTGCTA
60S ribosomal protein subunit
F: AGACCAGACAAAGCCAGAAGAC R: GCTGAAGTAAGTCCACGCATT
40S ribosomal protein subunit
F: GGCTACACCCACAGTTCTACCT R: CCAACGGTCTGGATAACAAA
E%: PCR amplifying efficiency; R2: Pearson correlation coefficients of the standard curves
21
2
Tm (°C)
Amplicon length (bp)
E%
R
65
342
97.4
0.999
65
383
98.0
0.994
65
117
96.9
0.998
65
230
96.7
0.998
65
142
97.6
0.990
65
247
104.6
0.995
65
224
97.1
0.998
65
337
97.8
1.000
65
371
98.2
0.999
65
149
95.9
1.000
Experiment Developmental
Developmental stages
Specise
Tissue
Gene Name
Stability value
Anaylsis Method
Crassostrea gigas
Integrity organism
Glyceraldehyde-3-phosphate
1.034
geNorm
Table 2. Validated HKGs under given experimental condition in mollusc species.
22
from oocyte to D-veliger
process
Different stages of the reproductive cycle Different stages of gametogenesis Different levels of haemocyte ploidy status Response to
Exposure to harbor
stress
pollution Bacteria challenging (Vibrio splendidus) Parasite challenging (Bonamia ostreae) Virus challenging (Ostreid Herpesvirus)
Tissue
dehydrogenase Elongation factor 1 α
0.885
Ovary
Elongation factor 1 α
0.930
Testis
18s ribosomal RNA
0.430
Mytilus edulis
Male gonads
Elongation factor 1
0.306
Mya arenaria
Haemocyte
Elongation factor 1
1.064
Ribosomal protein S18
1.186
Gill
Elongation factor 1 α
0.113
Digestive gland
Ribosomal protein L7
0.110
Pecten maximus L.
Mytilus spp. Mya arenaria
Haemocyte
0.330
Ribosomal protein S18 Elongation factor 1 α
Haemocyte
Ostrea edulis
Crassostrea gigas
Larvae Mantle(spat)
Crassostrea gigas
Glyceraldehyde-3-phosphate
0.049
NormFinder geNorm
NormFinder
Ribosomal protein L7
0.823
Ribosomal protein S18
0.822
Elongation factor 1 α
0.6
Ribosomal protein L5
0.172
NormFinder
Glyceraldehyde-3-phosphate
0.213
dehydrogenase
0.041
NormFinder geNorm
adp-ribosylation factor 1
Table 3. Stability values of ten candidate HKGs under designed experimental conditions calculated by NormFinder Experimental condition Tissue
geNorm
geNorm
cyclophilin
Various tissues
Gene symbol
NormFinder
0.046
dehydrogenase
Haliotis rufescens
distribution
Elongation factor 1
BestKeeper
OsHV-1 challenge
Temperature
23
OsHV-1 challenge
distribution 18s Actin EF-1β S18 NADH GAPDH EF-1α RL7 EF-1γ RL15
2.859 2.833 1.817 1.649 0.612 0.423 0.658 0.520 0.583 0.425
stress M 1.001 0.214 0.409 0.353 0.382 0.271 0.151 0.241 0.154 0.180
HE 0.617 0.332 0.249 0.196 0.375 0.292 0.293 0.239 0.127 0.171
Gl 1.509 / 0.377 0.402 / / 0.502 0.825 0.774 0.537
DG 0.754 / 0.736 0.378 / / 0.798 0.895 0.193 1.042
AM 1.473 0.192 0.363 0.332 0.261 0.222 0.173 0.342 0.281 0.269
Ft 0.696 0.263 0.271 0.265 0.241 0.250 0.111 0.188 0.161 0.320
M 0.041 0.029 0.007 0.021 0.122 0.035 0.045 0.122 0.014 0.043
HE 0.075 0.405 0.137 0.064 0.148 0.383 0.180 0.119 0.254 0.026
at different temperature M 0.220 0.106 0.089 0.040 0.135 0.115 0.134 0.059 0.051 0.040
HE 0.129 0.590 0.128 0.154 0.145 0.385 0.183 0.075 0.066 0.138
Bold numbers showed the top three stable HKGs under each given conditions. “M”, “HE”, “GI”, “DG”, “AM” and “Ft” stand for mantle tissue, heamocyte, gill, digestive glands, adductor muscle and foot respectively.
24
Figure legends Fig.1. Heat map of candidate HKGs expressing level in haemocytes under experimental condition T1 and T2. The expression values of candidate HKGs were estimated by calculating the reads per kilo base per million mapped reads (RPKM) in transcriptome data of haemocytes, the RPKM was normalized to the lowest one among six samples for each candidate HKG. T1: OsHV-1 challenge at low temperature (12°C); T2: OsHV-1 challenge at high temperature (18°C). Fig.2. Amplicons of ten candidate HKGs in 2% agarose gel. M: DNA ladder marker; 1: Actin, 2: 18s, 3: EF-1α, 4: EF-1β, 5: EF-1γ, 6: GAPDH, 7: NADH, 8: RL7, 9: S18 and 10: RL15. Fig.3. Expression levels (Cq values) of ten candidate HKGs in various tissues. The tissues include haemocytes, mantle tissue, gills, digestive glands, adductor muscles and foots. Each box represents lower quartile/median/ upper quartile while the whiskers represent the maximum and minimum values from each data set. Fig.4. Expression stability of ten HKGs analyzed by geNorm among various tissues of ark clams. a: average expression stability values of ten HKGs ranking; a lower value indicates a more stable expression. b: determination of the optimal number of HKGs needed for normalization. The pairwise variation (Vn/n+1) was calculated using n and n+1 reference genes by geNorm software. A cut-off value for the pairwise variation was suggested at 0.15, below this value indicating an additional HKG was not required. According to the geNorm software manual, the proposed 25
0.15 value must not be taken as a too strict cut-off. Even more than three genes were used, the pairwise variation didn’t show significant decrease. Thus three best HKGs (RL7, EF-1a and EF-1γ) were suggested in this case. Fig.5. Expression stability of ten HKGs in ark clams analyzed by geNorm under OsHV-1 challenge. A, B, C, D, E and F show expression stability value ranking of ten HKGs under OsHV-1 challenge in ark clam’s haemocytes, mantle tissue, adductor muscles, foots, gills and digestive glands respectively. a, b, c, d, e and f correspondingly presents the results of pairwise variation analysis (Vn/n+1) in ark clam’s haemocytes, mantle tissue, gills, digestive glands, adductor muscles and foots respectively. In haemocytes, mantles, adductor muscles and foots, the values of pairwise variation V2/3 were all below the cut-off value (0.15), and the best pair HKGs were suggested for expression analysis. Five HKGs were suggested in gill with least pairwise variation value (V5/6), and three HKGs (V3/4) in digestive glands. Fig.6. Expression stability of ten HKGs in ark clams analyzed by geNorm under temperature stress or OsHV-1 challenge at different temperature. A and B show expression stability value ranking of ten HKGs under temperature stress in ark clam’s haemocytes and mantle tissue, C and D show that under OsHV-1 challenge at different temperature in ark clam’s haemocytes and mantle tissue respectively. a and b presents the results of pairwise variation analysis (Vn/n+1) in ark clam’s haemocytes and mantle tissue corresponding to A and B, c and d are corresponding to C and D respectively. The values of pairwise variation V2/3 were all below 0.15 26
in haemocytes and mantle tissue under above two conditions, indicating two HKGs could satisfy gene expression analyses. Fig.7. Expression stability values of ten HKGs calculated by NormFinder. Tissue, stability values among various tissues; haemocyte-TV and mantle-TV, stability values under OsHV-1 challenge at different temperature in haemocytes and mantle tissue; haemocyte-T and mantle-T, stability values under temperature stress in haemocytes and mantle tissue; Haemocyte-V, mantle-V, foot-V and adductor muscle-V, stability values under OsHV-1 challenge in ark clams’ haemocytes, mantle tissue, foots and adductor muscles respectively.
27
Fi g1
Fi g2 M
1
2
3
4
5
6
7
8
9
1 0
Fi g3
Fi g4
Fi g5
Fi g6
Fi g7
Supplemental Data
Fig.S1. Expression levels (Cq values) of ten candidate housekeeping genes in mixed tissue samples. The hollow bars illustrate the expression levels of the candidate housekeeping genes while the whiskers indicate standard deviations
28
ACT RL7
RL15
18s
S18
NADH
GAPDH
EF-1β
EF-1γ
EF-1α
29
Fig.S2. Melting curve analysis of ten candidate reference genes. The melting curves
were generated by the Bio-Rad Bio-Rad CFX MangerTM software (vision 2.1)
30
Experiment Developmental process
Developmental stages from oocyte to D-veliger Different stages of the reproductive cycle Different stages of gametogenesis Different levels of haemocyte ploidy status
Response to
Exposure to harbor
stress
pollution Bacteria challenging (Vibrio splendidus) Parasite challenging (Bonamia ostreae) Virus challenging (Ostreid Herpesvirus)
Tissue distribution
Specise
Forward primer (5’-3’)
Reverse primer (5’-3’)
Glyceraldehyde-3-phosphate Crassostrea gigas
Pecten maximus L.
Mytilus edulis
Mya arenaria
Mytilus spp. Mya arenaria
Ostrea edulis
(°C)
(bp)
56
208
Elongation factor 1 α
F:CGAGAAGGAAGCTGCTGAGATGG
R:ACAGTCAGCCTGTGAAGTTC
60
172
Elongation factor 1 α
F:AGGGCTCCTTCAAGTATGCCTG
R:TGAGCGGTCTCGAACTTCCAC
60
100
18s ribosomal RNA
F: TCGATGGTACGTGATATGCC
R: CGTTTCTCATGCTCCCTCTC
60
84
Elongation factor 1
F:CACCACGAGTCTCTCCCAGA
R:GCTGTCACCACAGACCATTCC
60
106
Elongation factor 1
F:GGCCTAGGTGTTTTCCATGA
R:GGTGGCTGTTGGTGTCATC
60
158
Ribosomal protein S18
F:GCCGGTTGTCTTTGTATGCT
R:AAGATTCCCGACTGGTTCCT
60
189
Elongation factor 1 α
F:ACCCAAGGGAGCCAAAAGTT
R:TGTCAACGATACCAGCATCC
60
212
Ribosomal protein L7
F: CAGAGACAGGCCAAGAAAGG
R: TGGGTAGCCCCATGTAATGT
60
227
Elongation factor 1
F:GGTGGCTGTTGGTGTCATC
R:GGCCTAGGTGTTTTCCATGA
60
158
Ribosomal protein S18
F:GCCGGTTGTCTTTGTATGCT
R: AAGATTCCCGACTGGTTCCT
60
189
F:GTCGCTCACAGAAGCTGTACC
R:CCAGGGTGGTTCAAGATGAT
60
162
F:TCCCGCTAGCATTCCTTG
R:TTGGCGCCTCCTTTCATA
60
108
F:TCCCAAGCCAAGGAAGGTTATGC
R:CAAAGCGTCCAAGGTGTTTCTCAA
60
-
F:GCCATCAAGGGTATCGGTAGAC
R:CTGCCTGTTAAGGAACCAGTCAG
60
-
F:TCACCAACAAGGACATCATTTGTC
R:CAGGAGGAGTCCAGTGCAGTATG
60
152
F:CGCCAATCCTTGTTGCTT
R:TTGTCTTGCCCTCTTGC
-
-
F:AAGAACATTAGTTTCACTGT
R:CTATCTACCACAAAGATTAAC
-
-
Glyceraldehyde-3-phosphate Ribosomal protein L7 Ribosomal protein S18 Ribosomal protein L5 Glyceraldehyde-3-phosphate
Crassostrea gigas
Amplicons
R:CCTCTACCACCACGCCAATCCT
dehydrogenase Crassostrea gigas
Tm
F:AGAACATCGTCAGCAACGCATCC
dehydrogenase
Elongation factor 1 α
Haliotis rufescens
Various tissues
Gene Name
dehydrogenase ADP-ribosylation factor 1
31
Table S1. Primer information for validated housekeeping genes in moll
32
Graphical abstract
33
Highlights:
Previous studies of housekeeping genes selection in mollusc species were reviewed.
Ten candidate housekeeping genes were selected and assessed in ark shells under given conditions in ark shells.
RL7 and EF-1a were the most stable genes across tissues in ark shells.
RL15 appeared more stable than other housekeeping genes under various experimental conditions.
Commonly used housekeeping genes, 18s and ACT, were more greatly influenced by most given experimental conditions in ark shells.
34