Aquatic Toxicology 155 (2014) 129–141

Contents lists available at ScienceDirect

Aquatic Toxicology journal homepage: www.elsevier.com/locate/aquatox

Correlation of transcriptomic responses and metal bioaccumulation in Mytilus edulis L. reveals early indicators of stress Helen C. Poynton ∗ , William E. Robinson, Bonnie J. Blalock, Robyn E. Hannigan School for the Environment, University of Massachusetts Boston, 100 Morrissey Blvd., Boston, MA 02125, United States

a r t i c l e

i n f o

Article history: Received 31 May 2014 Received in revised form 23 June 2014 Accepted 24 June 2014 Available online 30 June 2014 Keywords: Toxicogenomics Tissue residue approach Mytilus edulis Lead Cadmium Quantitative RT-PCR

a b s t r a c t Marine biomonitoring programs in the U.S. and Europe have historically relied on monitoring tissue concentrations of bivalves to monitor contaminant levels and ecosystem health. By integrating ‘omic methods with these tissue residue approaches we can uncover mechanistic insight to link tissue concentrations to potential toxic effects. In an effort to identify novel biomarkers and better understand the molecular toxicology of metal bioaccumulation in bivalves, we exposed the blue mussel, Mytilus edulis L., to sub-lethal concentrations (0.54 ␮M) of cadmium, lead, and a Cd + Pb mixture. Metal concentrations were measured in gill tissues at 1, 2, and 4 weeks, and increased linearly over the 4 week duration. In addition, there was evidence that Pb interfered with Cd uptake in the mixture treatment. Using a 3025 sequence microarray for M. edulis, we performed transcriptomic analysis, identifying 57 differentially expressed sequences. Hierarchical clustering of these sequences successfully distinguished the different treatment groups demonstrating that the expression profiles were reproducible among the treatments. Enrichment analysis of gene ontology terms identified several biological processes that were perturbed by the treatments, including nucleoside phosphate biosynthetic processes, mRNA metabolic processes, and response to stress. To identify transcripts whose expression level correlated with metal bioaccumulation, we performed Pearson correlation analysis. Several transcripts correlated with gill metal concentrations including mt10, mt20, and contig 48, an unknown transcript containing a wsc domain. In addition, three transcripts directly involved in the unfolded protein response (UPR) were induced in the metal treatments at 2 weeks and were further up-regulated at 4 weeks. Overall, correlation of tissue concentrations and gene expression responses indicates that as mussels accumulate higher concentrations of metals, initial stress responses are mobilized to protect tissues. However, given the role of UPR in apoptosis, it serves as an early indicator of stress, which once overwhelmed will result in adverse physiological effects. © 2014 Elsevier B.V. All rights reserved.

1. Introduction Historically, national and international “mussel watch” programs including the NOAA National Status and Trends (NS&T) have

Abbreviations: AOP, adverse outcome pathway; BDL, below detection limits; DET, differentially expressed transcript; ER, endoplasmic reticulum; ERAD, ER-associated protein degradation; EST, expressed sequence tag; FDR, false discovery rate; GEO, gene expression omnibase; GO, gene ontology; HSP, heat shock protein; ICP-MS, inductively coupled plasma mass spectrometer; mt, metallothionein (gene); MT, metallothionein (protein); NS&T, National Status and Trends; POP, persistent organic pollutant; RDX, 1,3,5-trinitro-1,3,5-triazacyclohexane; ROS, reactive oxygen species; RT-qPCR, quantitative reverse transcription PCR; TCDD, 2,3,7,8-tetrachlorodibenzo-p-dioxin; TEQ, TCDD toxic equivalents; TNT, 2,4,6trinitrotoluene; TRA, tissue residue approach; UPR, unfolded protein response. ∗ Corresponding author. Tel.: +1 617 287 7323. E-mail address: [email protected] (H.C. Poynton). http://dx.doi.org/10.1016/j.aquatox.2014.06.015 0166-445X/© 2014 Elsevier B.V. All rights reserved.

monitored ecosystem health through the chemical analyses of contaminant levels in marine bivalve tissues (Goldberg, 1986; Cantillo, 1998). As sessile, ecologically-important coastal species, marine bivalves, in particular the blue mussel Mytilus edulis L., bioaccumulate a broad spectrum of organic and inorganic contaminants (Kimbrough et al., 2008). They are good surrogates for monitoring water quality because their soft tissues integrate exposure over time, even when water concentrations are close to detectable levels or when contaminant levels are temporally variable (Viarengo and Canesi, 1991; Blackmore and Wang, 2003). Several attempts have been made over the last 30 years to link contaminant body burdens with toxicity. While also referred to as the Critical Body Residue approach, recent work has adopted the term tissue residue approach (TRA), a term that encompasses both organic and inorganic contaminants, individual tissues and whole bodies (Meador et al., 2008). The theoretical basis for this approach assumes that the concentration of a contaminant at the site of toxic

130

H.C. Poynton et al. / Aquatic Toxicology 155 (2014) 129–141

action (i.e. the effective dose) is proportional to the concentration of that contaminant in the whole body and individual tissues. The closer the measurement of “dose” is in proximity to the target site, the more consistent responses between individuals or species will be, because major sources of variability such as bioavailability, uptake, toxicokinetics, and sequestration are bypassed (Meador et al., 2008). Therefore, by measuring either whole body or specific tissue contaminant concentrations, the TRA should provide a more accurate dose than measuring water or sediment contaminant concentrations, and should be more closely correlated with adverse effects (McCarty and Mackay, 1993; Meador et al., 2008). The TRA has proven to work very well for some classes of toxicants including chlorophenols and nonpolar organic compounds that elicit narcosis. In these cases, the assumption that whole body residues are proportional to the biologically active concentration at the target site holds up well (Meador et al., 2008). However, due to the variety and efficiency of mechanisms that have evolved in various species for metal uptake, detoxification and internal sequestration, tissue residues of metals may not correspond well to target site concentrations (Chapman, 1997; Rainbow, 2002; Luoma and Rainbow, 2005). Therefore, understanding the relationship between tissue concentrations and the biologically active concentration of metal (BAM) is not straightforward (Vijver et al., 2004; Adams et al., 2011). In addition to chemical analysis, many monitoring programs have recently incorporated biological measures of environmental health (i.e. biomarkers). A variety of biochemical, physiological and behavioral biomarkers has been developed for mussels and other bivalves in order to identify contaminant exposure, toxicological effects, and susceptibility to additional stressors (Decaprio, 1997; Shaw et al., 2011; Schettino et al., 2012). For example, the United Kingdom’s “Ecoman” monitoring program incorporates multiple biomarkers to prioritize sites for hierarchical risk assessment (Galloway et al., 2006). Over the past 10 years, ‘omic approaches have begun to provide additional opportunities for identifying novel biomarkers and unraveling modes of action in non-model organisms (Snape et al., 2004; Poynton et al., 2007; Poynton and Vulpe, 2009; Veldhoen et al., 2012). Transcriptomic, proteomic, and metabolomic biomarkers are likely to be among the first responses to contaminant stress, preceding that of physiological and behavioral biomarkers. Because of their enormous potential value to marine monitoring programs, there has been much interest in developing ‘omic resources and conducting transcriptomic studies on indicator organisms such as blue mussels. Venier et al. (2006) described a cDNA microarray for Mytilus galloprovincialis containing 1700 unique cDNAs, which provided transcriptome profiles capable of discriminating between different pollutants and between reference and impacted sites. Since then a number of additional studies have investigated transcriptomic responses in Mytilus spp. when exposed to okadaic acid (Manfrin et al., 2010), different mixtures of pollutants (Dondero et al., 2010, 2011; Canesi et al., 2011), heat or salinity stress (Lockwood et al., 2010; Lockwood and Somero, 2011; Mohamed et al., 2014), immunological stress (Venier et al., 2011; Philipp et al., 2012) and ocean acidification (Hüning et al., 2013). These studies have contributed mechanistic insight into the stress response of mussels and provided valuable biomarkers for use in monitoring programs. Some studies have even developed strategies for incorporating the complex data inherent in transcriptomic studies into monitoring databases. For example, Shaw et al. (2011) used gene expression profiling along with a suite of biomarkers within the Mussel Expert System (MES) to rate the level of stress experienced in mussels collected within the Tamar River and Estuary. Used in concert with TRA, ‘omics methods should provide mechanistic information to link tissue concentrations to toxic effects. For example, because gene expression responds to some portion of the biologically active internal dose of a chemical, when used

Cd, Pb, Mixture Exposures 1, 2, and 4 week me points

Differenal gene expression analysis (57 seqs)

Gill ssue metal concentraons

Pearson correlaon analysis

Differenally expressed genes correlated with ssue residue Fig. 1. Overview of research approach.

in combination with TRA, transcriptomics may be able to relate tissue concentrations to the biologically active concentration. This may be particularly helpful for the study of metal accumulation where presently, a lack of mechanistic detail related to the sequestration and detoxification of metals has prevented researchers from precisely estimating biologically active metals (BAMs) from total metal concentrations (Rainbow, 2002; Adams et al., 2011). Through expression profiles of cellular toxicity, ‘omics’ methods may identify biochemical pathways that are adversely impacted by BAMs and also identify pathways that trigger metal sequestration and detoxification pathways (e.g. metallothionein and antioxidant production) that mitigates toxicity and/or builds up the nonbioavailable metal pool. In addition, gene expression is more likely to be correlated with tissue metal concentrations than with water or sediment exposure concentrations. Despite the potential advantages of combining transcriptomics and TRA, only a few studies to date have made correlations between tissue concentrations and gene expression responses (Nakayama et al., 2006; Gong et al., 2012). The objective of the present study was to identify differentially expressed genes, through a transcriptomic approach, that correlated with metal accumulation in M. edulis gill tissue. Our approach, outlined in Fig. 1, involved 1-, 2- and 4-week exposures to cadmium, lead, and a Cd + Pb mixture, and subsequent measurement of gill gene expression responses and the corresponding gill metal concentrations in individual mussels. We first identified genes that were differentially expressed in each of the treatments and timepoints. From this reduced gene set, we then performed Pearson correlation analysis to identify genes whose expression level correlated with metal concentrations in the gills of the same mussels. Our goal was to provide monitoring programs with gene expression biomarkers that respond in a dose-dependent manner to specific metal bioaccumulation. In addition, we wanted to uncover insights into the mode of action and cellular toxicity of these metals through transcriptomic analysis, to help improve upon the TRA for metals in bivalves. 2. Materials and methods 2.1. Experimental animals One hundred and twenty M. edulis between 50 and 60 mm in length were collected in January 2011, from Bourne, MA (41.7402EN, 70.6157EW). Historical data from this NS&T site

H.C. Poynton et al. / Aquatic Toxicology 155 (2014) 129–141

(Kimbrough et al., 2008) indicate that the site is less contaminated by metals and organic toxicants than other coastal areas in Massachusetts. Mussels were collected by hand, at approximately the +0.15 m tidal height (4.1 ◦ C, 35 PSU salinity), and transported in a cooler back to the University of Massachusetts Boston. Once in the lab, mussels were cleaned of epifauna/flora, rinsed in cold tapwater, and placed in a 69 l (18 gal) aquarium containing 5 ␮m filtered seawater at 11 ◦ C for acclimation (8 d). Mussels were fed Isochrysis galbana (tropical strain) and Thalassiosira weissflogii every two days during the acclimation period. Over the course of the acclimation period (and experiment), mussels were subjected to a 10 h:14 h light:dark cycle, and salinities ranging from 30 to 32 PSU. 2.2. Experimental treatment design At T0 , seven mussels were randomly selected for initial analysis, dissected, and gill tissue was preserved as described below. Ninety two of the remaining mussels were randomly selected, and divided into each of four 19 l treatment aquaria (with polypropylene liner to minimize metal adsorption) filled with 16 l of 5 ␮m-filtered seawater and maintained in an environmental room (11 ◦ C; 10 h light/14 h dark). Mussels were exposed to 0.54 ␮M of the chloride salt of each metal alone or in combination: cadmium (Sigma Aldrich, St. Louis, MO) (61.26 ␮g/l Cd), lead (Sigma Aldrich, St. Louis, MO) (111.68 ␮g/l Pb) and a Cd + Pb mixture (61.26 ␮g/l Cd + 111.68 ␮g/l Pb) or an unexposed control. While slightly above the nanomolar concentrations found in estuaries (Bruland et al., 1979; Boutier et al., 1993; Jiann et al., 2005), it is lower than the concentrations used in other recent laboratory bivalve exposures (Vercauteren and Blust, 1999; Lemoine et al., 2000; Ciocan and Rotchell, 2004; Liu et al., 2012). MINTEQA2 modeling (Allison et al., 1991) of the experimental parameters (11 ◦ C, 31 PSU salinity, DIC entered as 2.4 × 10−3 M carbonate and pH either fixed at 8.0 or allowed to equilibrate) indicated that the relative speciation of Cd, Pb and all other major (≥1%) components of seawater were no different in the four treatment aquaria. We assumed independent uptake (and toxicity) of each metal (with no metal interactions). Aquarium water was changed twice per week (8 times during the experiment) and temperature, salinity, and dissolved oxygen were monitored in each tank (Table S1). Mortalities were mixture > Pb = control (Tamhane’s T2, p ≤ 0.05) indicating that less Cd was accumulated by mussels exposed to the mixture of Cd and Pb, than the mussels that were exposed to Cd alone. Similar to the Cd bioaccumulation, gill Pb concentration also increased linearly (if considering weeks 1, 2, and 4) or slightly curvilinearly (if including T0 ) over time in the Pb and mixture treatments (repeated measures, p < 0.01; Table 1), although no significant accumulation was noted in the Cd or control groups. The pattern of mean Pb concentration in the treatment groups was slightly different from the Cd analysis: Pb = mixture > Cd = control. In this case, similar amounts of Pb were accumulated in mussel gills regardless of whether exposure was to Pb alone or the mixture of Pb and Cd. Of the six other metals that were measured in gill tissue (Ca, Cr, Cu, Fe, Ni and Zn), all of them showed slight increases over time

H.C. Poynton et al. / Aquatic Toxicology 155 (2014) 129–141

133

Table 1 Summarized results of repeated measures analysis, conducted on log2 -transformed data. Metal

Cd Pb Ca Cr Cu Fe Ni Zn

Significant differences (p < 0.05) Among sampling days

Day × exposure interaction

Among exposure groups

Yes Yes Yes Yes Yes Yes Yes Yes

No Yes No No No No No No

Yes Yes No No No No (p = 0.065) No No

(Table S4), likely due to slightly elevated concentrations of metals in Boston Harbor seawater compared to the Cape Cod Canal seawater at the mussel collection site. None of these mean metal concentrations exhibited any differences among the four exposure groups (repeated measures, p < 0.05; Table 1). By the end of the experiment,

Differences among sampling days

Tamhane’s T2 differences among treatments Control

Pb

Mix

Cd

7 < 14 < 28 7 < 14 < 28 7 < 14 = 28 7 < 14 = 28 7 < 14 = 28 7 < 14 = 28 7 < 14 = 28 7 < 14 < 28

A A A A A A A A

A B A A A A A A

B B A A A A A A

C A A A A A A A

these metals exhibited a 1.6 to 4.3 fold increase compared to the T0 mussels, and showed a clear sign of plateauing out at either week 2 or 4. In contrast, our metal spiked treatments exhibited a considerably higher increase by the last sampling day (i.e. a 149× increase in Cd and a 200× increase in Pb). There was also an indication that both Cd and Pb were taken up from the Boston Harbor seawater (based on a 1.8-fold increase in Cd concentrations in the control and Pbexposed treatments, and a 4.3-fold increase in Pb concentrations in the control and Cd-exposed treatments). This increase was not statistically significant (Table 1). However, when comparing metal uptake in just these two sets of groups, roughly 2.5× more Pb was taken up than Cd, a finding consistent with what was observed in our metal spiked treatments.

3.2. Transcriptomic analysis with M. edulis microarray To identify genes differentially expressed in response to Cd or Pb bioaccumulation, we constructed a custom oligonucleotide microarray using the available sequences in GenBank (December 2010), which consisted of 3250 unique sequences following assembly. A total of 57 sequences were differentially expressed by at least one of the treatments and time-points (see Table S5 for a complete list of differentially expressed genes and mean expression levels and Table S6 for expression levels for individual mussels). The number of differentially expressed genes for each treatment at 1 week, 2 weeks, and 4 weeks are shown in Venn diagrams (Fig. 3). Week 1 showed the largest amount of overlap with five sequences differentially expressed by all three treatments, but only 1 or 2 sequences unique to each treatment. At two weeks, the largest number of sequences were differentially expressed and the majority of the sequences were specific to one treatment. By week 4, there were only eight sequences in total that were significantly differentially expressed. Because the largest number of sequences were differentially expressed at 2 weeks, we performed hierarchical clustering on the week 2 data set. Clustering of samples (Fig. 4) revealed that all the mussels within each treatment group clustered tightly together, demonstrating that despite potential genetic and life history

Fig. 2. Bioconcentration of Cd (A) and Pb (B) in gill tissue of Mytilus edulis after 7, 14 and 28 days of exposure to 0.54 ␮M concentrations of each metal in seawater ( denotes Control mussels; 䊉 denotes Cd exposure; , Pb exposure; and+Cd + Pb mixture). Error bars are S.E.

week 1: 20 sequences

week 2: 30 sequences

Pb

Pb

Cd 2

4 1

8

1

Cd

1

8 3

0

8

2 0

week 4: 10 sequences

Cd

3 1 0 2 1 0

Pb

3

6

11

Mix

Mix

Mix

Fig. 3. Number of differentially expressed sequences in each treatment at 1, 2, and 4 weeks represented in Venn Diagrams. The number of sequences which were differentially expressed in multiple treatments is shown in the overlapping areas of the circles.

134

H.C. Poynton et al. / Aquatic Toxicology 155 (2014) 129–141

Log2 expression level 4.6

9.0

16.4

Mix6 Mix2 Mix5 Mix4 Cd5 Cd3 Cd2 Cd6 ctrl7 ctrl3 ctrl1 ctrl2 Pb1 Pb5 Pb4 Pb3

Seq ID Seq. Description GO:0016071 mRNA metabolic processes, GO:0006412 Translation 40s ribosomal protein s4 Contig387 small nuclear ribonucleoprotein f-like Contig414 protein mago nashi homolog AM881703.1 dna-directed rna polymerase ii subunit rpb7 AM880894.1 60s ribosomal protein l13a DQ206347.1 60s ribosomal protein l6 Contig416 ribosomal protein l35a AM877668.1 AM878981.1 ribosomal protein s10 GO:0006950 Response to Stress protein transport protein sec61 GT157760.1 GO:19101293 Nucleocide phosphate biosynthetic process atp synthase subunit mitochondrial-like Contig430 AM877854.1 atp synthase subunit, mitochondrial ribonucleotide reductase m2 polypeptide AM882410.1 GO:0046872: Metal ion binding, Metallothioneins metallothionein 10 Contig20 metallothionein 20 AJ005456.1 Other Metabolic Processes foot protein 1 Contig21 AM878505.1 peptidyl-prolyl cis-trans isomerase fkbp2-like Contig419 protein rft1 homolog senesc ence-ass ociated protein AM878224.1 unknown function selenoprotein m-like protein AM879173.1 no homology found AM878081.1 AM879586.1 no homology found no homology found AM879936.1 no homology found AM880205.1 no homology found AM882350.1 no homology found GT157850.1 no homology found AM882376.1 no homology found Contig55 no homology found AM879452.1 no homology found AM882397.1 no homology found Contig428

Fig. 4. Hierarchical clustering of individual mussels after two week exposure to Cd, Pb, mixture of Cd and Pb, or the no exposure control. The 34 sequences which were differentially expressed at 2 weeks are shown, grouped according to shared gene ontology (GO) terms. Log2 normalized expression values for all differentially expressed are represented by the color scale.

differences between the field collected mussels, they responded similarly when challenged with Cd, Pb, or a mixture of the metals. A similar separation of treatment groups is shown at 1 week and 4 weeks (Fig. S1). At week 2, the mixture exposure of Cd and Pb falls on a main branch with the Cd exposure, suggesting that at this time-point Cd has a large influence on the gene expression pattern in the mixture exposure. However, by 4 weeks, the mixture exposure is separated from all other treatments, suggesting that the combined effects of Pb and Cd are impacting gene expression levels differently. To better understand the potential toxicity of the metal treatments, the M. edulis sequences were annotated using Blast2GO. Blast2GO assists in the annotation of genes by mapping their functions to a consistent set of vocabulary, or Gene Ontology (GO) terms that describe the biological processes, cellular component,

or the molecular function of the gene products (http://www. geneontology.org/). Grouping genes according to their GO terms can help to provide an overview of which physiological or biochemical processes are perturbed by a particular treatment. Of the 3250 unique sequences, 752 sequences were annotated and mapped to GO terms. This included 26 (or 46%) of the differentially expressed sequences (Table S7). Differentially expressed sequences in both Table S5 and Fig. 4 were organized according to similar GO terms. This revealed that several transcripts involved in mRNA metabolic processes, response to stress, nucleotide phosphate biosynthetic processes, and metal ion binding (two metallothionein transcripts) were differentially expressed. These four biological processes were also found to be overrepresented in our data set through enrichment analysis (Table 2), suggesting that these are the major biological processes that were perturbed by the metal

H.C. Poynton et al. / Aquatic Toxicology 155 (2014) 129–141

135

Table 2 Overrepresented gene ontology terms within the differentially expressed gene list. Enrichment analysis was performed using Fisher’s Test for enrichment with Blast2GO. GO-ID

Term

Molecular function Peroxidase activity GO:0004601 GO:0016209 Antioxidant activity GO:0020037 Heme binding Hydrolase activity, acting on glycosyl GO:0016798 bonds ATPase activity, coupled to GO:0044769 transmembrane movement of ions, rotational mechanism Proton-transporting ATPase activity, GO:0046961 rotational mechanism GO:0019829 Cation-transporting ATPase activity ribonucleoside-diphosphate reductase GO:0004748 activity, thioredoxin disulfide as acceptor Biological processes Response to stress GO:0006950 GO:0006979 GO:0009132 GO:1901293 GO:0008380 GO:0006397 GO:0016071

Response to oxidative stress Nucleoside diphosphate metabolic process Nucleoside phosphate biosynthetic process RNA splicing mRNA processing mRNA metabolic process

p-value

No. of test

No. of Ref.

Test seqs.

0.001 0.001 0.037 0.037

4 4 4 3

9 10 19 19

AM879017, Contig48, Contig87, AM880265

0.016

2

4

0.003

2

1

0.036 0.035

2 1

7 0

0.024

8

101

0.001 0.001

4 2

10 0

AM879017, Contig48, Contig220, Contig87, GT157760, AM880894, DQ201827, AM880265 AM879017, Contig48, Contig87, AM880265 AM882410, Contig430

0.028

3

17

AM877854, AM882410, Contig430

0.007 0.010 0.040

3 3 4

9 11 35

AM881703, AM880894, Contig414

treatments. Sequences mapped to nucleoside phosphate biosynthetic process (GO:1901293) included ATP synthase genes (Contig 430 and AM877854), which were down-regulated in the mixture at 4 weeks, and ribonucleotide reductase (AM882410) which was upregulated at 2 weeks in all treatments. mRNA metabolic processes (GO:0016071) represented by a number of ribosomal proteins were up-regulated in the different treatments at 2 weeks. We were particularly interested in the genes that mapped the biological process, response to stress (GO:0006950; Table 2). Of these eight genes, four of them also mapped to response to oxidative stress (GO:0006979), with homology to wsc domaincontaining protein. Despite their homology to the same gene, the four sequences showed two distinct gene expression patterns. Two sequences were down-regulated by Cd and Pb treatment at 1 week, while the other two sequences were down-regulated in the mixture and Pb at 4 weeks. Another response to stress sequence, X-box binding protein (DQ201827), generally showed higher expression in Pb treatment, but this up-regulation was only significant at 4 weeks. Two other sequences (Contig 220 and GT157760) involved in stress response were up-regulated in all treatments. In summary, four major biological processes were overrepresented in the differentially expressed sequence list, indicating that they were affected by the metal treatments. The inclusion of GO term, response to stress, suggests that the metal exposure were inducing a cellular stress response, which will be discussed in more detail below.

AM879017, Contig48, Contig87 AM877854, Contig430

AM882410

AM881703, Contig387, AM880894, Contig414

3.3. Targeted gene expression analysis by RT-qPCR A subset of differentially expressed sequences were chosen for confirmatory RT-qPCR. Six sequences were chosen to represent genes that were specific to either Cd or Pb, which could potentially serve as biomarkers for each of the metals. In addition, contig48 was chosen for its specific down-regulation by the mixture exposure and sec61 (GT157760) as a stress response gene (GO:0006950) that responded across all three treatments. Results from the eight RT-qPCR analyses are shown in Table 3 and a comparison of RTqPCR and microarray results is provided in Fig. S2. Overall, the gene expression patterns are consistent between the RT-qPCR and microarray analysis, despite representing different mussels. All of the sequences identified as differentially expressed in the microarray analysis were shown to be significant by RT-qPCR, with the exception of AM879173, selenoprotein (p = 0.068 for Pb, week 2). Based on the RT-qPCR results, two Cd specific transcripts were identified as members of the metallothionein gene family, mt20 (AJ005456) and mt10 (contig20). mt20 was highly up-regulated by Cd and the mixture at both the 2 week and 4 week time-points, while mt10 was up-regulated by Cd and mixture at 2 weeks, but in the Cd exposure, the expression level decreased at 4 weeks. Several genes that were significantly differentially expressed by the Pb treatment were confirmed by RT-qPCR including FK506 binding protein, AM878505, and X-box binding protein, DQ201827.

Table 3 Mean gene expression ratios of 2 week and 4 week exposed mussels vs. initial expression levels. Gene expression was assayed using RT-qPCR for each of the genes listed (See Table S2 for primer sequences and efficiencies). Relative gene expression was determined using the standard curve method and normalized to 18S rRNA. Normalized gene expression values were divided by the initial expression levels to derive the following expression ratios. Ratios that were significantly different from control are noted (* p < 0.05). Seq. Name

Seq. ID

Week 2 Ctrl

MT-20 Unknown FK506 Selenoprotein Unknown MT-10 X-box binding Sec61

AJ005456 Contig 48 AM878505 AM879173 AM882376 Contig 20 DQ201827 GT157760

2.45 0.10 0.36 0.39 0.95 1.53 1.01 0.60

Week 4 Cd

Pb *

379 1.39 0.68 0.66 1.03 10.5* 1.23 1.18

1.32 0.19 1.13* 0.73 2.39* 3.03 1.60* 1.16*

Mix *

170 0.12 0.63 0.55 0.76 9.10* 1.62* 1.05*

Ctrl 1.42 4.25 0.50 0.54 0.69 1.29 1.08 1.04

Cd *

173 4.34 0.63 0.44 0.91 2.40 1.10 1.17

Pb

Mix

0.65 0.29 1.22 1.25 1.87 2.73 2.67 1.50

238* 0.15* 0.92 0.88 0.85 5.37* 1.29 1.88

136

H.C. Poynton et al. / Aquatic Toxicology 155 (2014) 129–141

Table 4 Summary table of Cd and Pb Pearson’s correlation coefficients (with significance p-values in parentheses). Only coefficients that have p ≤ 0.05 are listed. Cd, mix, Pb and control are the 4 exposure treatments. n = 16 for each treatment and includes T0 points. Log2 DEG

Log2 [Pb] Cd

AJ005456.1 AM878224.1 AM878298.1 AM878436.1 AM878827.1 AM879452.1 AM880085.1 AM880265.1 AM881552.1 AM882350.1 AM882376.1 AM882397.1 Contig 416 Contig20 Contig23 Contig428 Contig48 Log2 DEG

AJ005456.1 AM878436.1 AM878849.1 AM878981.1 AM879452.1 AM880085.1 AM880265.1 AM881552.1 AM882376.1 AM882397.1 Contig20 Contig23 Contig48

Mix

Pb

Control

+0.731 (0.003) −0.513 (0.042) +0.504 (0.047) +0.501 (0.048) −0.717 (0.002)

+0.577 (0.019)

−0.524 (0.037) +0.619 (0.011) −0.503 (0.047) +0.639 (0.008) +0.508 (0.044) −0.658 (0.006) +0.601 (0.014)

−0.511 (0.043) +0.797 (0.000) −0.638 (0.008) −0.506 (0.045) −0.572 (0.021) Log2 [Cd] Cd

Mix

+0.884 (0.00)

+0.805 (0.00) +0.503 (0.047)

Pb

Control

+0.501 (0.048) +0.530 (0.035) +0.527 (0.036) +0.585 (0.017) +0.610 (0.012) −0.515 (0.041) +0.681 (0.004)

3.4. Correlation between gene expression and metal tissue concentrations From our list of differentially expressed sequences (total of 57), we aimed to identify genes whose expression pattern correlated with gill tissue metal concentrations. We performed a correlation analysis and determined the Pearson correlation coefficient (r) for each sequence by comparing its log2 expression level in the gill vs. the log2 metal concentration in the gill in each mussel. Table 4 provides a list of all the sequences with a significant linear relationship (p < 0.05) between the gene expression level and the level of either Cd or Pb measured in the gill tissue. For Pb, only one sequence (AM878827) showed a significant relationship between Pb and gene expression in the Pb treated mussels, although 13 genes showed a correlation between Pb and gene expression in mussels from the mixture exposure. For Cd, there were seven genes significantly correlated with Cd concentrations in the Cd treated mussels, and all but one of these sequences were correlated with Cd concentration in the mixture treatment as well. We also performed linear regression analysis and Pearson Correlation Analysis on the log2 transformed gene expression data vs. the log2 transformed metal accumulation data for the mussels that were analyzed by RT-qPCR (see Fig. 5). The RT-qPCR analyses differed from the microarray analysis in three important ways. It (1) utilized three different mussels in each exposure compared with the four mussels used in the microarray analysis, (2) did not include the week 1 time-point, and (3) included all treatments in a single analysis, as opposed to separating out each exposure treatment and analyzing each separately. Despite these differences, three sequences with significant correlations in the microarray analysis also showed a significant correlations with the RT-qPCR data

−0.529 (0.035) +0.562 (0.023) −0.520 (0.039) +0.673 (0.004) −0.626 (0.009) +0.593 (0.016) +0.761 (0.001) −0.632 (0.009) −0.633 (0.008)

+0.575 (0.020)

including AJ005456 (mt20), contig 20 (mt10), and contig 48. Two sequences with no significant correlation in microarray analysis showed a significant correlation with RT-qPCR gene expression data (DQ201827, X-box binding protein; GT157760, Sec61). One sequence (AM882376) with a significant correlation between Cd concentrations and gene expression in the microarray analysis, showed a slight, but not significant correlation when applied to the RT-qPCR data (r = -0.33, p = 0.09). In summary, the correlation analysis identified several transcripts whose expression correlated with the bioaccumulation of metals in the mussel gill tissue. Although the correlation analysis did not always agree between the microarray and RT-qPCR data sets, three transcripts, AJ005456 (mt20), contig 20 (mt10), and contig 48, were found to be robustly correlated and provide potential biomarkers of metal bioaccumulation. 4. Discussion 4.1. Metal accumulation Our data suggest that Cd does not affect the accumulation of Pb, but that Pb reduces the bioconcentration of Cd in the gills of M. edulis. However, as summarized by Norwood et al. (2003) and Vijver et al. (2011) in their recent meta-analyses on metal uptake studies, previous reports have shown conflicting patterns of metal interaction. Four independent studies documented a reduction in Cd bioconcentration in mussels only at the highest concentrations of accompanying metal tested (6.1 ␮M Zn, Phillips, 1976; 7.6 ␮M Zn, Jackim et al., 1977; 3.6 ␮M Zn or 0.32 ␮M Cu, Elliott et al., 1986; ≥10 ␮M Zn, Vercauteren and Blust, 1999). Other studies saw no changes in bioaccumulation of Cd in the presence of other metals. Neither Pb, Zn, Cu, Hg nor Fe affected Cd

H.C. Poynton et al. / Aquatic Toxicology 155 (2014) 129–141

12

A. MT-20 log2 gene expression

log2 gene expresssion

4

r = 0.894**

2 0 -2 -4 -6 -8 -10

C. contig48 r = -0.441*

10 8 6 4 2 0 -2

0

2

4

6

8

10

12

0 6

8

B. MT-10

log2 gene expresssion

log2 gene expresssion

137

r = 0.552** 6

4

2

0

2

4

6

8

10

12

4

6

8

10

12

4

6

8

10

12

D. X-Box binding protein r = 0.482*

4

2 0

2

4

6

8

10

12

0

log2 [Cd ] (nmol/g dry wt) log2 gene expresssion

6

2

E. Sec 61 r = 0.406*

4

2 0

2

log2 [Pb ] (nmol/g dry wt) Fig. 5. Correlations between gene expression and metal accumulation determined by linear regression. Gene expression values (from RT-qPCR) were log2 transformed and compared with log2 transformed metal concentrations measured in mussel tissue. Each point represents one mussel. Pearson correlation coefficients are provided for each correlation. Only genes with significant correlations are shown (* p < 0.05; ** p < 0.005). (A) Metallothionein-20, AJ005456; (B) Metallothionein-10, Contig 20; (C) Contig 48; (D) X-box binding protein, DQ201827; (E) Sec61, GT157760.

uptake in isolated M. edulis gills (Carpene and George, 1981). In a study most similar to ours (Ritz et al., 1982), binary mixtures of Cd (0.44 ␮M) and Pb (0.48 ␮M) did not affect the uptake of either metal over a ten day period. In a 10-day four-metal mixture (Cd, Pb, Cu and Zn), bioconcentration was the same as in single metal treatments (Ritz et al., 1982). Our results, at a similar metal exposure concentration (0.54 ␮M), showed no significant differences in gill Cd and Pb concentrations among treatments after only 1-week exposure. This is not surprising given the degree of inter-animal variability and the small sample size (n = 7 mussels per treatment) when considering this single time point alone. Only by including the two additional exposure time points (2and 4-week) and by considering all the data using repeated measures analysis do we see the adverse impact that Pb had on Cd uptake.

4.2. Transcriptomics Because our study utilized field collected organisms, our first goal was to assess reproducibility of expression profiles across replicate animals exposed to the same treatments. As shown in Fig. 4 and Fig. S1, replicate mussels consistently cluster together, demonstrating that the toxicant exposures produce reproducible

gene expression changes even when field collected animals with unknown and possible diverse genetic backgrounds are utilized. The three time points when the transcriptome was assessed provided a different pattern of gene expression (Fig. 3). After 1 week of exposure, 15 sequences were differentially expressed across the three treatments, and a third of these transcripts were in common between the treatments, suggesting a similar physiological response. However, by 2 weeks, the gene expression profile shifted and 85% of the differentially expressed transcripts were unique to the individual treatments. Similarly, Wu and Wang (2010) studied metabolomic effects of Perna viridis exposed to Cu and Cd and found that at 1 week, the metabolic profiles overlapped between the individual treatments and a mixture, but they became more distinct by 2 weeks. Finally after 4 weeks of exposure, the fewest number of M. edulis sequences were differentially expressed and the majority of transcripts were specific to only one treatment. The pattern seen here, that with longer duration exposures the gene expression profile becomes more distinct, suggests that initially the cellular response to metal stress was similar. However, as the gill tissues bioaccumulated higher concentrations of the metals, and possibly incurred cellular damage due to the increasing concentrations, the response became more specific for each metal. Whether this pattern is limited to our study, or a common response to metal bioaccumulation needs to be determined with further studies.

138

H.C. Poynton et al. / Aquatic Toxicology 155 (2014) 129–141

4.2.1. Insights into modes of toxicity The concentrations chosen for the present study (0.54 ␮M) are far below those known to elicit acute mortality in mussel species (i.e. 24-h LC50 values for Cd at 14.2 ␮M and Pb at 21.7 ␮M in M. galloprovincialis; Vlahogianni and Valavanidis, 2007). However, sublethal effects have been shown at concentrations relevant to the present study including chronic effects to growth (Strömgren, 1982), decreased Ca2+ uptake, and gill filament injuries including inflammation, loss of cilia, and necrosis (Domouhtsidou and Dimitriadis, 2000; Sheir et al., 2013). Sublethal effects of Cd and Pb to gill tissue are thought to be mediated through lipid peroxidation and oxidative damage (Vlahogianni and Valavanidis, 2007). Transcriptomic analysis revealed that several genes involved in stress responses were differentially expressed, including transcripts with a role in the unfolded protein response (UPR) (Cao and Kaufman, 2012). These include three transcripts with homology to Sec61 (GT157760), X-box binding protein (DQ201827), and FK506 binding protein (AM878505), which displayed a similar, significantly correlated gene expression pattern across the exposures (Table S8). While expression of these transcripts was similar to control values at week 1, the transcripts were induced at week 2, and further up-regulated at 4 weeks (see Table S5, Fig. S2, and Fig. 5D and E). UPR provides a coordinated response to misfolded proteins in the endoplasmic reticulum (ER) leading to enhanced ER function, a reduction in overall protein synthesis, and the induction of ER–associated protein degradation (ERAD) (Travers et al., 2000; Cao and Kaufman, 2012). The X-box binding protein is a transcription factor which plays a critical role in ERAD through transcriptional regulation of multiple genes involved with ER function and secretatory pathways (Boyce and Yuan, 2006). Sec61 is a component of the translocon involved in transport of newly synthesized polypeptides into the ER. It is up-regulated during UPR (Travers et al., 2000) and may also be involved in the ERAD by trafficking damaged ER-associated proteins out of the ER (Boyce and Yuan, 2006; Hampton and Sommer, 2012). The third transcript involved in UPR and induced by our metal exposures was FK506 binding protein, a peptidyl-prolyl cis–trans isomerase. This protein is localized to the ER lumen where it is involved in protein folding and trafficking, and is induced under conditions that lead to protein mis-folding (Bush et al., 1994). UPR can be induced by a variety of cellular stresses including Ca2+ depletion and oxidative stress, both of which have been shown to occur upon Cd exposure in marine bivalves (Chora et al., 2008; Sheir et al., 2013), and in mammalian renal tubular cells, Cd has been shown to directly induce UPR (Kitamura and Hiramatsu, 2010). Finally, Contig 220 was upregulated in response to the mixture exposure. This transcript codes for ubiquinin conjugating enzyme E2 that targets proteins for degradation. Chora et al. (2008) showed that clams exposed to 0.36 ␮M Cd have increased protein ubiquination, and later confirmed that Cd caused a decrease in total number of proteins in gills (Chora et al., 2009). These findings are consistent with an activation of the UPR, which includes increased protein degradation of misfolded proteins and decrease in protein translation (Cao and Kaufman, 2012). Effects were also seen in transcripts related to mRNA metabolism and translation. This class of transcripts was not coordinated and displayed fluctuating patterns of gene expression across the treatments. For example, several transcripts coding ribosomal proteins were upregulation by Pb, but other transcripts involved in mRNA metabolism had a mixed response by Cd and the mixture (see Table S5). Other studies have shown that heavy metal exposure affects mRNA metabolism related transcripts in M. galloprovincialis suggesting that disruption to mRNA metabolism may be a general response to heavy metal exposure (Venier et al., 2006; Dondero et al., 2011).

Finally, four transcripts coding for wsc domain containing proteins were differentially expressed (see Table S5). These genes are best characterized in yeast as cell wall stress sensors (Lodder et al., 1999) and therefore homology of four transcripts to this class of proteins resulted in an overrepresentation of the GO term, response to oxidative stress (GO:0006979). However, the role of these proteins in metazoans is less clear. The wsc domain may be involved in carbohydrate binding (pfam01822), but it is found in many proteins of diverse function. Therefore, subscribing it to a specific function in oxidative stress is not possible. 4.2.2. Expression profiles of metal mixtures Because Cd is a more potent toxicant compared with Pb to mussels (Yap et al., 2004; Vlahogianni and Valavanidis, 2007), and we exposed the mussels to equimolar concentrations of the metals, we expected the gene expression profiles from the mixture exposed mussels to be dominated by the influence of Cd and resemble the Cd expression profile. This is clearly the case after 2 weeks, where hierarchical clustering of gene expression profiles separates the mixture and Cd exposure on a distinct branch (Fig. 4 and Fig. S1). However, this pattern is not seen at the other time-points. In fact, after 4 weeks of exposure, the mixture is separated from the other treatments (Fig. S1). It is interesting to note that at 4 weeks, the mussels exposed to the mixture have accumulated only about half the concentration of Cd compared with the mussels from the Cd only treatment (Fig. 2A), and this difference is significant (Table 1). The differences in Cd bioaccumulation at 4 weeks may partially account for the divergence of gene expression profiles. In two binary mixture studies in mussels, the gene expression profiles of the mixtures following short-term exposures (96-h), clustered with one of the two chemicals, suggesting that similarly to our 2 week time-point, one chemical had a greater influence on overall gene expression (Canesi et al., 2011; Dondero et al., 2011). 4.3. Correlation of gene expression with metal bioaccumulation 4.3.1. Correlation analysis To identify transcripts that correlated with the bioaccumulation of Pb and Cd in gill mussel tissue, we employed a two step process outlined in Fig. 1. We first identified genes which were differentially expressed in at least one of the treatments and time-points. Using this sub-set of 57 transcripts, we conducted Pearson correlation analysis to identify those transcripts whose expression profile significantly correlated in a linear, dose-dependent fashion with either Pb or Cd concentrations (Table 4). While other patterns of dose–response may also be important (e.g. curvilinear or threshold dose response), we felt that a linear dose-dependent response would prove to be the most useful for identifying practical, yet robust, biomarkers. Our approach was similar to other recent studies targeting genes which were associated with contaminant concentrations in tissues. Studies on wild birds and wildlife used correlation analysis to identify genes associated with persistent organic pollutant (POP) accumulation and total 2,3,7,8-tetrachlorodibenzop-dioxin (TCDD) toxic equivalents (TEQs) (Nakayama et al., 2006; Hirakawa et al., 2011). In contrast to our study, these studies applied correlation analysis over their entire microarray data set. Gong et al. (2012) demonstrated a bioinformatic method for identifying gene biomarkers that corresponded to tissue residue concentrations of 2,4,6-trinitrotoluene (TNT) and 1,3,5-trinitro1,3,5-triazacyclohexane (RDX) in earthworms (Gong et al., 2012). Their approach involved an optional filtering step for differentially expressed genes, similar to the present study, and then used LASSO or random forest regression analysis to identify transcripts correlated with tissue residues. Although our approach was less

H.C. Poynton et al. / Aquatic Toxicology 155 (2014) 129–141

sophisticated, it provides a simple method for identifying transcripts correlated to tissue concentrations. In addition, the initial filtering step for differentially expressed genes reduces the total number of correlations that have to be calculated and thus the computational power needed to identify transcripts likely correlated with bioaccumulation. 4.3.2. Metallothionein genes Two metallothionein transcripts, metallothionein-10 (mt10) and metallothionein-20 (mt20), were highly correlated with Cd tissue concentrations. Metallothioneins (MTs) are low molecular weight, soluble proteins with a high cysteine content and the characteristic repeated Cys-X-Cys motif. The cysteine residues participate in the binding of many essential and nonessential metals in cysteine tetrathiolate clusters. It is thought that MTs are involved in metal sequestration and detoxification (Viarengo et al., 1999b), and possibly protect against oxidative stress (Viarengo et al., 1999a). In the present study, mt20 induction occurred at week 1 and its expression level in both the Cd and mixture treatments was highly induced (100–400 fold) throughout the entire exposure duration (Table S5, Fig. S2). In contrast, mt10 expression was only significantly different from the control treatments at week 2 in the Cd exposure and in week 2 and week 4 in the mixtures. The level of induction was also less impressive (3–10 fold) compared with mt20 at these time-points. These results are consistent with other studies that have shown distinct expression patterns for mt10 and mt20 (Barsyte et al., 1999; Dondero et al., 2005). For example, Banni et al. (2007) investigated patterns of mt10 and mt20 gene expression in M. galloprovincialis through in situ field exposures. They found that mt10 was constitutively expressed, and induced by a broad range of essential and nonessential metals while mt20 was expressed at very low levels under normal conditions and induced up to 1000-fold in the presence of Cd (Banni et al., 2007). Because of their inducibility upon metal exposure, there has been interest in incorporating MT as a biomarker of exposure in monitoring programs both at the protein (Viarengo et al., 1999b) and transcript levels (Lemoine and Laulier, 2003; Banni et al., 2007). 4.3.3. Identification of biomarkers of metal bioaccumulation One aim of our study was to identify transcripts that correlated with metal bioaccumulation in the gill of M. edulis that could serve as potential biomarkers of bioaccumulation. M. edulis is the primary mussel species used on the Northeast United States by the NOAA NS&T (Kimbrough et al., 2008) and the identification of novel biomarkers of bioaccumulation could help to prioritize sites for further investigation and provide mechanistic insight into the health of the surveyed mussels. Our approach successfully identified both major isoforms of mt as differentially expressed by Cd treatments and also correlated with Cd bioaccumulation in gill tissues. Because several studies have suggested the use of mt gene induction in field studies or monitoring programs (Banni et al., 2007; Franzellitti et al., 2010), their identification served as an important proof-ofprinciple for our approach. In addition to mt genes, we also identified novel transcripts that correlated with metal bioaccumulation. For example, contig 48 correlated with Pb bioaccumulation. Although its function is unknown, its significant correlation with Pb tissue concentrations makes this a potentially important biomarker of Pb contamination. Finally, we identified three transcripts with roles in the UPR that were differentially expressed in the metal treatments. These transcripts respond specifically to ER stress, but play a similar role to heat shock proteins (HSP) in responding to damaged proteins. Indeed, several HSPs are also involved in UPR signaling (Gorman et al., 2012). HSPs have also been proposed by many for use as biomarkers in monitoring programs; however, given the large number of both natural and unnatural stresses that induce them, attributing their

139

up-regulation to pollution is difficult (Feder and Hofmann, 1999; Fabbri et al., 2008). Perhaps, the UPR related transcripts identified in the present study, may provide a more specific indication of pollution induced stress; although, this will need to be verified in future studies. However, as many have demonstrated (Dondero et al., 2006; Viarengo et al., 2007; Franzellitti et al., 2010), using the identified transcripts together as a suite of biomarkers will likely provide a better prediction of metal bioaccumulation than a single transcript alone.

4.3.4. Integration of TRA and transcriptomic approaches The goal of using the TRA is to reduce variability in toxicity thresholds by approximating the concentration of a toxic chemical at its target site of action. However, total body concentrations or even specific tissue concentrations of contaminants do not accurately represent the biologically active concentration of many classes of chemicals, especially for metals where mechanisms for cellular sequestration, metabolism, of detoxification exist (Vijver et al., 2004; Adams et al., 2011). Therefore, a complementary method is needed which can provide mechanistic insight into when the biologically active concentration has reached a toxicity threshold and is causing adverse effects. By providing holistic signatures of the molecular physiology of cells (Poynton and Vulpe, 2009), transcriptomic approaches have the potential to provide the mechanistic information needed to determine when toxicity thresholds are breached by the biologically active chemical. In the present study, we found a correlation between tissue concentrations of both Pb and Cd and transcripts involved in the UPR. The coordinated up-regulation of three important genes in this response, indicated that in the mussels with higher accumulated metal concentrations, damage to proteins was occurring. However, the induction of these genes may indicate that compensation mechanisms were effectively utilized to reduce cellular toxicity. The high induction of mt20 in the Cd and mixture exposures also suggests that Cd was being sequestered, potentially reducing its biologically active concentration. Overall, these responses suggest that the mussels in the present study had mobilized the appropriate compensation mechanisms and are in within the “resistant” stage of the stress response, according to Seyle’s General Adaptation Syndrome Model (Selye, 1950). However, further accumulation of Cd and Pb could cause these mechanisms to become overwhelmed, resulting in toxicity. For example, while UPR serves to protect the cells from protein damage, UPR signaling instead leads to apoptosis during prolonged stress (Gorman et al., 2012). Indeed, in mammalian cells exposed to Cd, apoptosis mediated through the UPR has been shown to occur (Kitamura and Hiramatsu, 2010). Therefore, in combination, UPR and mt induction may serve as early indicators of metal induced stress in monitoring programs, although additional validation of this biomarker suite is needed. As genomic resources for environmentally important species grows, and annotation of these genomes increases, future studies will be able to more thoroughly characterize the physiological state of tissues through transcriptomics. The next steps will be to develop adverse outcome pathways (AOPs) to link molecular responses to adverse effects at the organismal and population levels (Ankley et al., 2010; Kramer et al., 2011; Hutchinson et al., 2013). Differential expression of transcripts mapped to an AOP will be a strong predictor of toxic effects providing evidence for when tissue concentrations have reached a toxicological threshold. As analytical chemical methods also improve, constructing a dose–response relationship of a chemical at its target site of action may soon be possible, which coached within an AOP will provide the basis for a truly mechanistic risk assessment model (Gutsell and Russell, 2013).

140

H.C. Poynton et al. / Aquatic Toxicology 155 (2014) 129–141

5. Supporting Information Additional supporting figures and tables are available online as supplemental files. These include water quality data for exposures (Table S1), complete sequences for ESTs and contigs used to construct the M. edulis microarray (Table S2), primer sequences for qPCR (Table S3), total metal concentrations in mussel gill tissues (Table S4), complete list of differentially expressed ESTs and expression values (Table S5 and S6), Blast2GO annotation of differentially expressed genes (Table S7), and Pearson correlation analysis from qPCR data (Table S8). Hierarchical clustering results for differentially expressed genes at each week (Fig. S1) and comparison of microarray and qPCR expression levels for select genes (Fig. S2) are also provided.

Acknowledgements Support for this study was provided by faculty start-up funds and a Healey Research Grant from the University of Massachusetts Boston to H.C.P. This publication is partially the result of research sponsored by The MIT Sea Grant College Program, under NOAA Grant number NA14OAR4170077, MIT SG Project Number 2014R/RCM-37 to H.C.P and W.E.R.

Appendix A. Supplementary data Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/ j.aquatox.2014.06.015.

References Adams, W.J., Blust, R., Borgmann, U., Brix, K.V., DeForest, D.K., Green, A.S., Meyer, J.S., McGeer, J.C., Paquin, P.R., Rainbow, P.S., Wood, C.M., 2011. Utility of tissue residues for predicting effects of metals on aquatic organisms. Integr. Environ. Assess. Manage. 7, 75–98. Allison, J.D., Brown, D.S., Novo-Gradac, K.J., 1991. MINTEQA2/PRODEFA2, a geochemical assessment model for environmental systems. In: U.S. EPA/600/3-91/021, Athens, GA. Ankley, G.T., Bennett, R.S., Erickson, R.J., Hoff, D.J., Hornung, M.W., Johnson, R.D., Mount, D.R., Nichols, J.W., Russom, C.L., Schmieder, P.K., 2010. Adverse outcome pathways: a conceptual framework to support ecotoxicology research and risk assessment. Environ. Toxicol. Chem. 29, 730–741. Banni, M., Dondero, F., Jebali, J., Guerbej, H., Boussetta, H., Viarengo, A., 2007. Assessment of heavy metal contamination using real-time PCR analysis of mussel metallothionein mt10 and mt20 expression: a validation along the Tunisian coast. Biomarkers 12, 369–383. Barsyte, D., White, K.N., Lovejoy, D.A., 1999. Cloning and characterization of metallothionein cDNAs in the mussel Mytilus edulis L. digestive gland. Comp. Biochem. Phys. C 122, 287–296. Benjamini, Y., Hochberg, Y., 1995. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J. R. Stat. Soc. B Met. 57, 289–300. Bevington, P.R., Robinson, D.K., 1992. Data Reduction and Error Analysis for the Physical Sciences, 2nd ed. WCB/McGraw Hill, Boston, MA. Blackmore, G., Wang, W.X., 2003. Comparison of metal accumulation in mussels at different local and global scales. Environ. Toxicol. Chem. 22, 388–395. Boutier, B., Chiffoleau, J., Auger, D., Truquet, I., 1993. Influence of the Loire River on dissolved lead and cadmium concentrations in coastal waters of Brittany. Estuarine, Coastal Shelf Sci. 36, 133–145. Boyce, M., Yuan, J., 2006. Cellular response to endoplasmic reticulum stress: a matter of life or death. Cell Death Differ. 13, 363–373. Bruland, K.W., Franks, R.P., Knauer, G.A., Martin, J.H., 1979. Sampling and analytical methods for the determination of copper, cadmium, zinc, and nickel at the nanogram per liter level in sea water. Anal. Chim. Acta 105, 233–245. Bush, K., Hendrickson, B., Nigam, S., 1994. Induction of the FK506-binding protein, FKBP13, under conditions which misfold proteins in the endoplasmic reticulum. Biochem. J. 303, 705–708. Canesi, L., Negri, A., Barmo, C., Banni, M., Gallo, G., Viarengo, A., Dondero, F., 2011. The organophosphate chlorpyrifos interferes with the responses to 17␤-estradiol in the digestive gland of the marine mussel Mytilus galloprovincialis. PLoS One 6, e19803. Cao, S.S., Kaufman, R.J., 2012. Unfolded protein response. Curr. Biol. 22, R622–R626. Carpene, E., George, S.G., 1981. Absorption of cadmium by gills of Mytilus edulis (L.). Mol. Physiol. 1, 23–34.

Cantillo, A., 1998. Comparison of results of mussel watch programs of the United States and France with worldwide mussel watch studies. Mar. Pollut. Bull. 36, 712–717. Chapman, P.M., 1997. Is bioaccumulation useful for predicting impacts. Mar. Pollut. Bull. 34, 282–283. Chora, S., McDonagh, B., Sheehan, D., Starita-Geribaldi, M., Roméo, M., Bebianno, M.J., 2008. Ubiquitination and carbonylation as markers of oxidative-stress in Ruditapes decussatus. Mar. Environ. Res. 66, 95–97. Chora, S., Starita-Geribaldi, M., Guigonis, J.-M., Samson, M., Roméo, M., Bebianno, M.J., 2009. Effect of cadmium in the clam Ruditapes decussatus assessed by proteomic analysis. Aquat. Toxicol. 94, 300–308. Ciocan, C.M., Rotchell, J.M., 2004. Cadmium induction of metallothionein isoforms in juvenile and adult mussel (Mytilus edulis). Environ. Sci. Technol. 38, 1073–1078. Conesa, A., Gotz, S., Garcia-Gomez, J.M., Terol, J., Talon, M., Robles, M., 2005. Blast2GO: a universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics 21, 3674–3676. Decaprio, A.P., 1997. Biomarkers: coming of age for environmental health and risk assessment. Environ. Sci. Technol. 31, 1837–1848. Domouhtsidou, G.P., Dimitriadis, V.K., 2000. Ultrastructural localization of heavy metals (Hg, Ag, Pb, and Cu) in gills and digestive gland of mussels, Mytilus galloprovincialis (L.). Arch. Environ. Contam. Toxicol. 38, 472–478. Dondero, F., Banni, M., Negri, A., Boatti, L., Dagnino, A., Viarengo, A., 2011. Interactions of a pesticide/heavy metal mixture in marine bivalves: a transcriptomic assessment. BMC Genomics 12, 195. Dondero, F., Dagnino, A., Jonsson, H., Capri, F., Gastaldi, L., Viarengo, A., 2006. Assessing the occurrence of a stress syndrome in mussels (Mytilus edulis) using a combined biomarker/gene expression approach. Aquat. Toxicol. 78 (Suppl. 1), S13–S24. Dondero, F., Negri, A., Boatti, L., Marsano, F., Mignone, F., Viarengo, A., 2010. Transcriptomic and proteomic effects of a neonicotinoid insecticide mixture in the marine mussel (Mytilus galloprovincialis, Lam.). Sci. Total Environ. 408, 3775–3786. Dondero, F., Piacentini, L., Banni, M., Rebelo, M., Burlando, B., Viarengo, A., 2005. Quantitative PCR analysis of two molluscan metallothionein genes unveils differential expression and regulation. Gene 345, 259–270. Elliott, N., Swain, R., Ritz, D., 1986. Metal interaction during accumulation by the mussel Mytilus edulis planulatus. Mar. Biol. 93, 395–399. Fabbri, E., Valbonesi, P., Franzellitti, S., 2008. HSP expression in bivalves. Invertebr. Surviv. J. 5, 135–161. Feder, M.E., Hofmann, G.E., 1999. Heat-shock proteins, molecular chaperones, and the stress response: evolutionary and ecological physiology. Annu. Rev. Physiol. 61, 243–282. Franzellitti, S., Buratti, S., Donnini, F., Fabbri, E., 2010. Exposure of mussels to a polluted environment: insights into the stress syndrome development. Comp. Biochem. Phys. C 152, 24–33. Galloway, T.S., Brown, R.J., Browne, M.A., Dissanayake, A., Lowe, D., Depledge, M.H., Jones, M.B., 2006. The ECOMAN project: a novel approach to defining sustainable ecosystem function. Mar. Pollut. Bull. 53, 186–194. Goldberg, E.D., 1986. The mussel watch concept. Environ. Monit. Assess. 7, 91–103. Gong, P., Loh, P.R., Barker, N.D., Tucker, G., Wang, N., Zhang, C., Escalon, B.L., Berger, B., Perkins, E.J., 2012. Building quantitative prediction models for tissue residue of two explosives compounds in earthworms from microarray gene expression data. Environ. Sci. Technol. 46, 19–26. Gorman, A.M., Healy, S.J.M., Jaeger, R., Samali, A., 2012. Stress management at the ER: regulators of ER stress-induced apoptosis. Pharmacol. Ther. 134, 306–316. Gutsell, S., Russell, P., 2013. The role of chemistry in developing understanding of adverse outcome pathways and their application in risk assessment. Toxicol. Res. 2, 299–307. Hampton, R.Y., Sommer, T., 2012. Finding the will and the way of ERAD substrate retrotranslocation. Curr. Opin. Cell Biol. 24, 460–466. Hirakawa, S., Imaeda, D., Nakayama, K., Udaka, M., Kim, E.Y., Kunisue, T., Ogawa, M., Matsuda, T., Matsui, S., Petrov, E.A., Batoev, V.B., Tanabe, S., Iwata, H., 2011. Integrative assessment of potential effects of dioxins and related compounds in wild Baikal seals (Pusa sibirica): application of microarray and biochemical analyses. Aquat. Toxicol. 105, 89–99. Hüning, A.K., Melzner, F., Thomsen, J., Gutowska, M.A., Krämer, L., Frickenhaus, S., Rosenstiel, P., Pörtner, H.-O., Philipp, E.E., Lucassen, M., 2013. Impacts of seawater acidification on mantle gene expression patterns of the Baltic Sea blue mussel: implications for shell formation and energy metabolism. Mar. Biol. 160, 1845–1861. Hutchinson, T.H., Lyons, B.P., Thain, J.E., Law, R.J., 2013. Evaluating legacy contaminants and emerging chemicals in marine environments using adverse outcome pathways and biological effects-directed analysis. Mar. Pollut. Bull. 74, 517–525. Jackim, E., Morrison, G., Steele, R., 1977. Effects of environmental factors on radiocadmium uptake by four species of marine bivalves. Mar. Biol. 40, 303–308. Jiann, K.-T., Wen, L.-S., Santschi, P.H., 2005. Trace metal (Cd, Cu, Ni and Pb) partitioning, affinities and removal in the Danshuei River estuary, a macro-tidal, temporally anoxic estuary in Taiwan. Mar. Chem. 96, 293–313. Kimbrough, K.L., Johnson, W.E., Lauenstein, G.G., Cristensen, J.D., Apeti, D.A., 2008. An assessment of two decades of contaminant monitoring in the nation’s coastal zone. In: NOAA Technical Memorandum NOS NCCOS 74, Silver Spring, MD.

H.C. Poynton et al. / Aquatic Toxicology 155 (2014) 129–141 Kitamura, M., Hiramatsu, N., 2010. The oxidative stress: endoplasmic reticulum stress axis in cadmium toxicity. Biometals 23, 941–950. Kramer, V.J., Etterson, M.A., Hecker, M., Murphy, C.A., Roesijadi, G., Spade, D.J., Spromberg, J.A., Wang, M., Ankley, G.T., 2011. Adverse outcome pathways and ecological risk assessment: bridging to population-level effects. Environ. Toxicol. Chem. 30, 64–76. Lemoine, S., Bigot, Y., Sellos, D., Cosson, R., Laulier, M., 2000. Metallothionein isoforms in Mytilus edulis (Mollusca, Bivalvia): complementary DNA characterization and quantification of expression in different organs after exposure to cadmium, zinc, and copper. Mar. Biotechnol. 2, 195–203. Lemoine, S., Laulier, M., 2003. Potential use of the levels of the mRNA of a specific metallothionein isoform (MT-20) in mussel (Mytilus edulis) as a biomarker of cadmium contamination. Mar. Pollut. Bull. 46, 1450–1455. Liu, F., Wang, D.Z., Wang, W.X., 2012. Cadmium-induced changes in trace element bioaccumulation and proteomics perspective in four marine bivalves. Environ. Toxicol. Chem. 31, 1292–1300. Lockwood, B.L., Sanders, J.G., Somero, G.N., 2010. Transcriptomic responses to heat stress in invasive and native blue mussels (genus Mytilus): molecular correlates of invasive success. J. Exp. Biol. 213, 3548–3558. Lockwood, B.L., Somero, G.N., 2011. Transcriptomic responses to salinity stress in invasive and native blue mussels (genus Mytilus). Mol. Ecol. 20, 517–529. Lodder, A.L., Lee, T.K., Ballester, R., 1999. Characterization of the Wsc1 protein, a putative receptor in the stress response of Saccharomyces cerevisiae. Genetics 152, 1487–1499. Luoma, S.N., Rainbow, P.S., 2005. Why is metal bioaccumulation so variable? Biodynamics as a unifying concept. Environ. Sci. Technol. 39, 1921–1931. Manfrin, C., Dreos, R., Battistella, S., Beran, A., Gerdol, M., Varotto, L., Lanfranchi, G., Venier, P., Pallavicini, A., 2010. Mediterranean mussel gene expression profile induced by okadaic acid exposure. Environ. Sci. Technol. 44, 8276–8283. McCarty, L.S., Mackay, D., 1993. Enhancing ecotoxicological modeling and assessment. Body residues and modes of toxic action. Environ. Sci. Technol. 27, 1718–1728. Meador, J.P., McCarty, L.S., Escher, B.I., Adams, W.J., 2008. 10th anniversary critical review: the tissue-residue approach for toxicity assessment: concepts, issues, application, and recommendations. J. Environ. Monit. 10, 1486–1498. Mohamed, B., Hajer, A., Susanna, S., Caterina, O., Flavio, M., Hamadi, B., Aldo, V., 2014. Transcriptomic responses to heat stress and nickel in the mussel Mytilus galloprovincialis. Aquat. Toxicol. 48, 104–112. Nakayama, K., Iwata, H., Kim, E.Y., Tashiro, K., Tanabe, S., 2006. Gene expression profiling in common cormorant liver with an oligo array: assessing the potential toxic effects of environmental contaminants. Environ. Sci. Technol. 40, 1076–1083. Norwood, W., Borgmann, U., Dixon, D., Wallace, A., 2003. Effects of metal mixtures on aquatic biota: a review of observations and methods. Hum. Ecol. Risk Assess. 9, 795–811. Philipp, E.E., Kraemer, L., Melzner, F., Poustka, A.J., Thieme, S., Findeisen, U., Schreiber, S., Rosenstiel, P., 2012. Massively parallel RNA sequencing identifies a complex immune gene repertoire in the lophotrochozoan Mytilus edulis. PLoS One 7, e33091. Phillips, D., 1976. The common mussel Mytilus edulis as an indicator of pollution by zinc, cadmium, lead and copper. I. Effects of environmental variables on uptake of metals. Mar. Biol. 38, 59–69. Poynton, H.C., Varshavsky, J.R., Chang, B., Cavigiolio, G., Chan, S., Holman, P.S., Loguinov, A.V., Bauer, D.J., Komachi, K., Theil, E.C., Perkins, E.J., Hughes, O., Vulpe, C.D., 2007. Daphnia magna ecotoxicogenomics provides mechanistic insights into metal toxicity. Environ. Sci. Technol. 41, 1044–1050. Poynton, H.C., Vulpe, C., 2009. Ecotoxicogenomics: emerging technologies for emerging contaminants. J. Am. Water Resour. Assoc. 45, 83–96. Rainbow, P.S., 2002. Trace metal concentrations in aquatic invertebrates: why and so what? Environ. Pollut. 120, 497–507. Ritz, D., Swain, R., Elliott, N., 1982. Use of the mussel Mytilus edulis planulatus (Lamarck) in monitoring heavy metal levels in seawater. Mar. Freshwater Res. 33, 491–506.

141

Rozen, S., Skaletsky, H., 2000. Primer3 on the WWW for general users and for biologist programmers. Methods Mol. Biol. 132, 365–386. Schettino, T., Caricato, R., Calisi, A., Giordano, M.E., Lionetto, M.G., 2012. Biomarker approach in marine monitoring and assessment: new insights and perspectives. Open Environ. Sci. 6, 20–27. Selye, H., 1950. Stress and the general adaptation syndrome. Br. Med. J. 1, 1383. Shaw, J.P., Dondero, F., Moore, M.N., Negri, A., Dagnino, A., Readman, J.W., Lowe, D.R., Frickers, P.E., Beesley, A., Thain, J.E., Viarengo, A., 2011. Integration of biochemical, histochemical and toxicogenomic indices for the assessment of health status of mussels from the Tamar Estuary, UK. Mar. Environ. Res. 72, 13–24. Sheir, S.K., Handy, R.D., Henry, T.B., 2013. Effect of pollution history on immunological responses and organ histology in the marine mussel Mytilus edulis exposed to cadmium. Arch. Environ. Contam. Toxicol. 64, 701–716. Snape, J.R., Maund, S.J., Pickford, D.B., Hutchinson, T.H., 2004. Ecotoxicogenomics: the challenge of integrating genomics into aquatic and terrestrial ecotoxicology. Aquat. Toxicol. 67, 143–154. Strömgren, T., 1982. Effect of heavy metals (Zn, Hg, Cu, Cd, Pb, Ni) on the length growth of Mytilus edulis. Mar. Biol. 72, 69–72. Travers, K.J., Patil, C.K., Wodicka, L., Lockhart, D.J., Weissman, J.S., Walter, P., 2000. Functional and genomic analyses reveal an essential coordination between the unfolded protein response and ER-associated degradation. Cell 101, 249–258. Veldhoen, N., Ikonomou, M.G., Helbing, C.C., 2012. Molecular profiling of marine fauna: integration of omics with environmental assessment of the world’s oceans. Ecotoxicol. Environ. Saf. 76, 23–38. Venier, P., De Pitta, C., Pallavicini, A., Marsano, F., Varotto, L., Romualdi, C., Dondero, F., Viarengo, A., Lanfranchi, G., 2006. Development of mussel mRNA profiling: can gene expression trends reveal coastal water pollution? Mutat. Res. 602, 121–134. Venier, P., Varotto, L., Rosani, U., Millino, C., Celegato, B., Bernante, F., Lanfranchi, G., Novoa, B., Roch, P., Figueras, A., Pallavicini, A., 2011. Insights into the innate immunity of the Mediterranean mussel Mytilus galloprovincialis. BMC Genomics 12, 69. Vercauteren, K., Blust, R., 1999. Uptake of cadmium and zinc by the mussel Mytilus edulis and inhibition by calcium channel and metabolic blockers. Mar. Biol. 135, 615–626. Viarengo, A., Burlando, B., Cavaletto, M., Marchi, B., Ponzano, E., Blasco, J., 1999a. Role of metallothionein against oxidative stress in the mussel Mytilus galloprovincialis. Am. J. Physiol. 277, R1612–R1619. Viarengo, A., Burlando, B., Dondero, F., Marro, A., Fabbri, R., 1999b. Metallothionein as a tool in biomonitoring programmes. Biomarkers 4, 455–466. Viarengo, A., Canesi, L., 1991. Mussels as biological indicators of pollution. Aquaculture 94, 225–243. Viarengo, A., Lowe, D., Bolognesi, C., Fabbri, E., Koehler, A., 2007. The use of biomarkers in biomonitoring: a 2-tier approach assessing the level of pollutantinduced stress syndrome in sentinel organisms. Comp. Biochem. Phys. C 146, 281–300. Vijver, M.G., Elliott, E.G., Peijnenburg, W.J., De Snoo, G.R., 2011. Response predictions for organisms water-exposed to metal mixtures: a meta-analysis. Environ. Toxicol. Chem. 30, 1482–1487. Vijver, M.G., Van Gestel, C.A., Lanno, R.P., Van Straalen, N.M., Peijnenburg, W.J., 2004. Internal metal sequestration and its ecotoxicological relevance: a review. Environ. Sci. Technol. 38, 4705–4712. Vlahogianni, T., Valavanidis, A., 2007. Heavy-metal effects on lipid peroxidation and antioxidant defence enzymes in mussels Mytilus galloprovincialis. Chem. Ecol. 23, 361–371. Wu, H., Wang, W.X., 2010. NMR-based metabolomic studies on the toxicological effects of cadmium and copper on green mussels Perna viridis. Aquat. Toxicol. 100, 339–345. Yap, C.K., Ismail, A., Omar, H., Tan, S.G., 2004. Toxicities and tolerances of Cd, Cu, Pb and Zn in a primary producer Isochrysis galbana and in a primary consumer Perna viridis. Environ. Int. 29, 1097–1104.

Correlation of transcriptomic responses and metal bioaccumulation in Mytilus edulis L. reveals early indicators of stress.

Marine biomonitoring programs in the U.S. and Europe have historically relied on monitoring tissue concentrations of bivalves to monitor contaminant l...
1MB Sizes 3 Downloads 4 Views