RNA-Seq analysis of the toxicant-induced transcriptome of the marine diatom, Ceratoneis closterium

MARGEN-00179; No of Pages 9 Marine Genomics xxx (2013) xxx–xxx

Contents lists available at ScienceDirect

Marine Genomics journal homepage: www.elsevier.com/locate/margen

4 5 6 7

Sharon E. Hook a,⁎, Hannah L. Osborn a, Francesca Gissi a, Philippe Moncuquet b, Natalie A. Twine c, Marc R. Wilkins c, Merrin S. Adams a a b c

CSIRO Land and Water, Locked Bag 2007, Kirrawee, NSW 2232 Australia CSIRO Computational Informatics, Acton, ACT, Australia NSW Systems Biology Initiative, School of Biotechnology and Biomolecular Sciences, University of New South Wales, Sydney, NSW 2052, Australia

8

i n f o

a b s t r a c t

P

a r t i c l e

Article history: Received 27 September 2013 Received in revised form 29 November 2013 Accepted 18 December 2013 Available online xxxx

D

Diatoms are of enormous ecological importance as they account for as much as 20% of global primary production, yet they are still understudied from a genomic perspective. The benthic diatom Ceratoneis closterium is wellcharacterized from an ecotoxicological perspective including its use in ecotoxicological risk assessments and investigating the mode-of-action of metal toxicity. However, this organism has little sequence information available. In this study, 454 pyrosequencing of the stressor-responsive transcriptome was undertaken. These transcripts could be used to characterize general physiological processes such as photosynthesis and respiration, as well as to enable a description of the ecotoxicogenomic responses of this organism. After a 96 h exposure to the concentration of toxicant that inhibited growth rate by 10% (IC10) for the following common coastal contaminants: ammonia, copper, crude oil and simazine (a photosystem II inhibiting herbicide), diatom cells were harvested for RNA extraction and their transcriptomes characterized via 454 pyrosequencing. This resulted in 1.25 million reads, which were assembled into 4768 contigs, when contigs encoding rRNA were removed. More than 80% of the remaining contigs had an ortholog in the BLASTx protein databases. These contigs represented 1660 unique transcripts. The role of these transcripts in stress response, as well as photosynthesis and respiration is discussed. Overall, this study greatly enhances the genomic information available for this important taxonomic group. © 2013 Published by Elsevier B.V.

C

T

E

Keywords: Diatom Ecotoxicology Photosynthesis Next generation sequencing Gene expression Stress response

R

E

9 10 11 12 13 14 16 15 17 18 19 20 21 22 23 24

F

3Q1

O

2

RNA-Seq analysis of the toxicant-induced transcriptome of the marine diatom, Ceratoneis closterium

R O

1

44

R

43

1. Introduction

46

The field of ecotoxicogenomics, which has emerged over the last decade, has led to many advances in ecotoxicology. In particular, ecotoxicogenomics has enabled the detection of responses to contaminants at lower concentrations and over shorter time periods than is possible in traditional tests (Poynton and Vulpe, 2009). It has also allowed the exploration of mode of action of toxicants and interaction effects, enabling better prediction of how contaminants will affect organisms in “real world” scenarios (Villeneuve and Garcia-Reyero, 2011). However, even with the advent of next generation sequencing technology (reviewed in Mehinto et al., 2013), the majority of transcriptomic studies are still performed with either microarrays or qPCR, both of which require a priori knowledge of the organisms' genome. Marine diatoms account for as much as 20% of global primary production and make large contributions to oceanic biogeochemical cycling

49 50 51 52 53 54 55 Q5 56 57 58 59 60

U

47 48

N C O

45

⁎ Corresponding author. Tel.: +61 2 9710 6839. E-mail address: [email protected] (S.E. Hook).

25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 42 41

(Falkowski et al., 2004; Parker et al., 2008). Despite their ecological importance, there have been comparatively few ecotoxicogenomic studies performed on diatoms (e.g. Hook and Osborn, 2012; Osborn and Hook, 2013) in part because of the paucity of sequence information available. At the time of writing (June 2013), there were only approximately 32,000 mRNA sequences in NCBI's GenBank nucleotide database. Most of these sequences belonged to only three species: Thalassiosira pseudonana (11,713), Phaeodactylum tricornutum (10,480), and Thalassiosira oceanica (9614). The remainder of species had less than 100 sequences associated with them in the nucleotide database. In Australia, the temperate marine diatom Ceratoneis closterium has been utilised to assess the toxicity and bioavailability of contaminants in coastal environments for over 20 years. This species is also often referred to by its previous names, Nitzschia closterium and Cylindrotheca closterium. Toxicity tests with C. closterium measure chronic toxicity to the microalga as inhibition in growth rate over a 72-h exposure to contaminants (Stauber et al., 1994). Due to the algals' sensitivity to a range of contaminants, the diatom has been used in a range of ecotoxicological studies including toxicity assessments of sewage effluent (Adams et al., 2008), landfill leachates (Binet et al., 2001), sediments (MorenoGarrido et al., 2003) and, investigation of the mechanisms of copper

1874-7787/$ – see front matter © 2013 Published by Elsevier B.V. http://dx.doi.org/10.1016/j.margen.2013.12.004

Please cite this article as: Hook, S.E., et al., RNA-Seq analysis of the toxicant-induced transcriptome of the marine diatom, Ceratoneis closterium, Mar. Genomics (2013), http://dx.doi.org/10.1016/j.margen.2013.12.004

61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 Q6

2.1. Algae cultures

112

123

The marine diatom Ceratoneis closterium (Ehrenb.) W. Smith was originally obtained from CSIRO Collection of Living Microalgae, Marine and Atmospheric Research, Hobart, Australia. This species was recently reinstated as Ceratoneis closterium, but had been previously been known as N. closterium, and is maintained in the culture collection under its previous name (Ehrenb.) Reimann and J.C. Lewin (http:// www.csiro.au/Organisation-Structure/National-Facilities/AustralianNational-Algae-Culture-Collection.aspx). Axenic cultures were maintained in silica containing f medium with the iron and trace elements halved (Guillard and Ryther, 1962), at 21 ± 2 °C, under a 12:12 h light:dark cycle, with an irradiance of 75 ± 5 μmol s−1 m−2 (TL 40 W cool white fluorescent lighting).

124

2.2. Chemical analysis

125

134 135

The water accommodated fraction (WAF) of crude oil was manufactured as outlined by Aurand and Coelho (2005), and used for toxicity tests as described in Osborn and Hook (2013). WAF samples were analysed using US EPA standard methods by the National Measurement Institute (Pymble, NSW, Australia). Analyte measurement and the methods used are described in Osborn and Hook (2013), and include dissolved polycyclic aromatic hydrocarbons (PAHs) and total petroleum hydrocarbons (TPHs). As the copper concentration used (2.5 μg/L) was close to the level of detection for our system (2 μg/L), we elected to work with nominal concentrations of copper. Nominal concentrations were also used for ammonia and simazine.

136

2.3. Toxicity tests

137

The chronic toxicity effects of copper, ammonia, crude oil (as WAF) and simazine on C. closterium were determined using 72 h growthrate inhibition bioassays following the procedure outlined by Stauber et al. (1994) and Franklin et al. (2005). Toxicity tests were carried out

89 90 91 92 93 94 95 96 Q8 97 98 99 100 101 102 103 104 105 106 107 108

F

111

88

O

2. Materials and methods

86 87

R O

110

84 85

in 250 mL borosilicate glass Erlenmeyer flasks (Bacto, NSW Australia) coated with a Coatasil silanizing solution (Coatasil, BDH) to prevent contaminants from adsorbing to the glass. Toxicity tests on the WAF were carried out in 20 mL non-silanized vials with Teflon-lined lids. Methanol was used as the carrier solvent in simazine tests, and a methanol control (0.05% w/v) was included in each toxicity test. Test solutions of copper and ammonia were prepared from stock solutions of copper sulphate and ammonium chloride respectively. The toxicity tests were carried out in filtered (0.45 μm) seawater (pH 8.2 ± 0.1, salinity 35 ± 2‰) collected from Cronulla, NSW, Australia, and supplemented with nitrate (15 mg NO−3 L−1) and phosphate (1.5 mg PO3 −4 L−1) to maintain exponential growth over the 72 h test duration. The volume of test solution in each flask was 50 mL except for tests conducted in vials where the test volume was 10 mL. Three replicates per treatment and a control were prepared. Cells in exponential growth phase (5–6 days old) were washed three times with filtered seawater by centrifugation (2500 rpm, 7 min, rotor radius 17 cm, Spintron, Melbourne, VIC, Australia) to remove nutrient rich media. Each test vessel was inoculated with 1 × 104 cells mL−1 of algal suspension. Tests were conducted under the same culture conditions described above with the light intensity increased to 107 ± 17 μmol s−1 m−2 (TL 40 W cool white fluorescent lighting). Sub-samples were taken after 24, 48 and 72 h to determine cell densities using the Becton Dickinson FACSCalibur flow cytometer, and growth rate was calculated as described in Franklin et al. (2005). Statistical endpoints for the concentrations required to reduce growth by 10% (IC10), and for the concentrations required to reduce growth by 50% (IC50) were determined using linear interpolation (Toxcalc. Scientific Tidepool Software Version 5.0.23, San Francisco, CA, USA).

E

109

toxicity to marine microalgae (Stauber and Forence, 1987; Levy et al., 2007). Toxicity tests with this species have also been adapted to identify the toxicants causing toxicity in waste waters (Hogan et al., 2005). We aim to apply modern molecular toxicological approaches to our work with the diatom Ceratoneis closterium to improve our understanding of the effect and mechanistic action of low concentrations of contaminants; however, this alga has no genomic information available. Next generation sequencing has been used to determine the sequence of coding genes as well as measure differences in gene expression in other organisms that lacked prior sequence information (e.g. Lowe et al., 2011; Zeng et al., 2011; Craft et al., 2010; Hudson, 2008; Burns et al., 2013). Here we describe the application of NGS to an ecotoxicogenomic study of the diatom Ceratoneis closterium. The aims were 1) to describe the transcriptome of the diatom Ceratoneis closterium, both to serve as a reference in future RNA-Seq projects, as well as to enable the construction of microarrays and qPCR assays and 2) to record changes in the transcriptome when the diatom is exposed to contaminants. The organism was cultured in control conditions, or with separate exposure to copper, ammonia, oil, and simazine. These contaminants were selected because they were regularly identified as contaminants of concern (Osborn and Hook, 2013) and have different modes of toxic action to C. closterium. The concentrations of copper, ammonia, oil and simazine selected were those that caused a 10% inhibition in algal growth rate, that is, the concentration at which a small effect on chronic toxicity is observed. The results of this study provide the only transcriptomic resource available for this ecologically and toxicologically significant species. They also provide a description of potential marker genes for future environmental studies.

P

82 Q7 83

S.E. Hook et al. / Marine Genomics xxx (2013) xxx–xxx

D

2

126 127 128 129 Q9 130 131 132 133

138 139 140

C

E

R

R

121 122

O

119 120

C

117 118

N

115 116

U

113 114

T

2.4. RNA extractions

141 142 143 144 Q10 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170

Algae were exposed to each stressor for 96 h at the concentrations shown in Table 1. The growth period was extended to 96 h so that a sufficient quantity of RNA could be extracted. Cells were collected as described in Osborn and Hook (2013) and centrifuged to produce a concentrated pellet. The overlying seawater was removed and the remaining pellet immediately immersed in liquid nitrogen, then stored at −80 °C. Five replicates per treatment were harvested. RNA extractions were carried out as described in Hook and Osborn (2012). Briefly, cells were lysed in a 4 M guanidine thiocyanate buffer (containing 0.75 M sodium citrate, 10% (v/v) sarcosyl and 0.36% (v/v) β mercaptoethanol), vortexed in a lysing matrix E (Molecular BioScience) then RNA was extracted using a combined TRIzol (Invitrogen) and RNeasy (Qiagen) procedure. Samples were then TURBO DNase treated (Applied Biosystems). RNA quantity and quality was verified with a bioanalyser (Agilent). Ribosomal RNA was removed using Invitrogen's Ribominus kit for plants.

171 172

2.5. Preparation of cDNA libraries and pyrosequencing

187

Construction of cDNA libraries and pyrosequencing was conducted at the Ramaciotti Centre for Gene Function Analysis, University of New South Wales, Kensington, NSW, Australia. Libraries were prepared as described in “cDNA Rapid Library Preparation Method Manual — GS FLX Titanium Series — October 2009 (Rev. Jan 2010)”, except for the following modifications: RNA was fragmented for 125 s; enzymatic reactions and size selection were performed on the SpriWorks (Beckman Coulter). Samples were multiplexed using MID-labelled primers (multiplex identifiers) RL5, RL6, RL7, RL8, and RL9. Combined cDNAs were then sequenced on a picotitre plate using the GS-FLX platform (454, Roche, Maryland, USA). Emulsion PCR (emPCR) titrations were carried out as described in the “emPCR Method Manual — Lib-L SV — GS FLX Titanium Series — October 2009 (Rev. Jan 2010).” Bulk emPCRs were carried out as described in the “emPCR Method Manual — Lib-L LV — GS FLX Titanium Series — October 2009 (Rev. Jan

188 189

Please cite this article as: Hook, S.E., et al., RNA-Seq analysis of the toxicant-induced transcriptome of the marine diatom, Ceratoneis closterium, Mar. Genomics (2013), http://dx.doi.org/10.1016/j.margen.2013.12.004

173 174 175 176 177 178 179 180 181 182 183 184 185 186

190 191 192 193 194 195 196 197 198 199 200 201 202

S.E. Hook et al. / Marine Genomics xxx (2013) xxx–xxx Table 1 Chronic toxicity of contaminants to C. closterium determined by growth inhibition.

t1:3

Toxicant

Units

N

IC10b

t1:4 t1:5 t1:6 t1:7

Copper Ammonia Simazine WAF

μg/L mg N/La μg/L g/L

4 4 4 4

2.2 4.9 310 1.0

t1:8 t1:9

a

7.1 ± 4.5 9.3 ± 1.9 N1000 10 ± 3

2.5 4 300 1

Total ammonia. The IC10 concentration inhibits growth rate by 10%, the IC50 concentration inhibits growth rate by 50% ± 1 standard deviation.

209

2.6. Data analysis

Prior to assembly, the sequencing reads were initially filtered for quality using the CSIRO instance of Galaxy (Goecks et al., 2010; 212 Giardine et al., 2005). The first 17 bp were trimmed from the 5′ end of 213 the sequence (to ensure the removal of MIDs and adaptors). After trim214 ming, a sliding window approach was used to filter out the reads with 215 low quality. The window size was 20, the minimum quality score was 216 20 (Blankenberg et al., 2010) and the number of bases to exclude was 217 1. Filtered reads were assembled using the MIRA algorithm (version 218 0.0.4) (Chevreux et al., 1999), as implemented in the CSIRO instance 219 of Galaxy. Reads were deposited into NCBI's databases with the acces220 sion number (to be determined). This Transcriptome Shotgun Assembly 221 project has been deposited at DDBJ/EMBL/GenBank under the accession 222 GAPN00000000. The version described in this paper is the first version, 223 Q11 GAPN01000000. The data are also available at: Hook et al. (in press): 224 raw reads and final assembled contigs for pennate diatom project. v1. 225 Q12 CSIRO. Data collection http://dx.doi.org/10.4225/08/525BCD6E6F145. The assembled contigs were submitted to Blast2Go v2.5.0 (http:// 226 227 www.blast2go.org/; Conesa et al., 2005) for sequence annotation and 228 Gene Ontology mapping. Contig sequences were compared to the 229 non-redundant protein database ‘nr’ using BLASTx (within the Blast2GO 230 tool). Sequences were assigned an annotation if the BLASTx search 231 resulted in a match with an e-value b 1 × 10−3 (however, the software 232 did not return any matches with e values greater than 9 × 10− 7). 233 Sequence annotation was made more parsimonious using the GO 234 SLIM annotations and the “plant” option (McCarthy et al., 2006). Taxo235 nomic distribution of results from a BLASTn search (version 0.0.11) 236 (conducted via the CSIRO instance of Galaxy) was obtained using the 237 “metagenomic analysis” tools within Galaxy. 238 Contigs were filtered to remove rRNA using CLC Genomics work239 bench (CLC Genomics workbench v6.0.1. CLC Bio A/S, Denmark). 240 Contigs were mapped to the 04_08_2012 downloads from SILVA 241 (Quast et al., 2013) using default parameters, and only “unmapped” 242 contigs (i.e. those that did not encode ribosomal subunits) were used 243 for differential expression analysis. Differential expression was calculat244 ed using CLC Genomics workbench. To calculate mapping percentages, 245 the “map reads to reference” tool was used with default parameters. 246 The results, given in reads per kilobyte mapped, were divided by the 247 number of reads for each treatment, to account for the different number 248 of reads created for each library. The “RNA-Seq Analysis” tool within CLC 249 Genomics was used with default parameters to map reads to the assem250 bly. This mapping was used to derive expression levels. Since only one 251 cDNA library was created per treatment, no statistical analysis could 252 be performed. Since read counts of less than 10 may not be a reliable 253 Q13 measure of true transcript abundance (Anders and Huber, 2010), 254 contigs with smaller read counts were not included in the differential 255 expression analysis.

3.1. Toxicity tests

257

P

D

E

T

C

E

R

R

N C O

256

The IC10 and IC50 values for exposure of C. closterium to copper, ammonia, simazine and WAF are presented in Table 1. These contaminants were chosen because they are commonly found in coastal ecosystems and known to cause toxicity via different modes of toxic action (Osborn and Hook, 2013). Simazine caused significant inhibition in growth, compared to controls, at concentrations greater than 500 μg/L. The IC50 concentration could not be calculated for simazine because it was greater than the highest concentration tested (1000 μg/L). The solvent control (methanol) showed no effect on growth inhibition (data not shown). WAFs of crude oil caused significant inhibition in growth at concentrations greater than 1 g/L. IC10 level exposures were used to maximise the toxicant-specific transcriptomic responses without eliciting a generalized stress response (e.g. Poynton et al., 2007; Hook et al., 2010).

210 211

U

3. Results and discussion

O

207 208

2010).” Sequencing was performed on a GS FLX (Roche) using the GS FLX Titanium Sequencing Kit XLR70 as described in the “Sequencing Method Manual — GS FLX Titanium Series — October 2009 (Rev. November 2010).” Image and signal processing were performed with the GS Sequencer, contained within Roche's 454 system sequencing software package version 2.6.

205 206

1.3 0.6 70 0.5

Exposure concentration for RNA extractions

R O

203 204

b

± ± ± ±

IC50b

F

t1:1 t1:2

3

258 259 260 261 262 263 264 265 266 267 268 269 270 271

3.2. Assembly statistics

272

Five different cDNA libraries for C. closterium were generated: a control library, and four separate libraries for exposure to ammonia, copper, oil (WAF) and simazine. These samples were multiplexed with specialized primers and run on one Roche 454 pyrosequencing plate. The 454 sequencing of the C. closterium pooled transcriptome generated 1,248,932 raw reads. The assembly data is summarized in Table 2. After trimming to remove low quality bases from reads MIRA — as implemented in the CSIRO instance of galaxy — assembled only approximately 20% of the available reads (263,589 reads of the 1.22 million reads that passed QA) into 12,779 contigs with an average length of 906 bp. Although rRNA was removed before cDNA was generated for sequencing, evidently this process was inefficient as more than 50% of our contigs encoded rRNA. Removing contigs that mapped to rRNA (as represented by the SILVA database) (Quast et al., 2013) resulted in 4768 non-rRNA contigs remaining with an average length of 1146 bp. The length distribution for this subset of contigs is shown in Fig. 1.

273 274

3.3. Sequence annotation

290

When the rRNA was removed by mapping contigs to the SILVA database (Quast et al., 2013), only 4769 remained suggesting that much of the fragmented assembly encoded ribosomal subunits. However, 3997 of these contigs had a protein ortholog in the BLAST database and 3131 have a nucleotide ortholog. From the nucleotide ortholog output, 1660 of the nucleotide hits are unique. However, as shown in Fig. 2 not all of the contigs align to an ortholog in BLAST for the entire length of the sequence. In addition, there were very few unique orthologs, as shown in Table 3. Furthermore, as these contigs aligned to less than 2000 different coding products, it suggests that the final assembly is still fragmented. Of the unique BLAST nucleotide orthologs, 1599 had taxonomic information associated with them. As shown in Fig. 3, most sequences aligned best to eukaryotes (1515 — 95%). Within the eukaryotes, most

291

Please cite this article as: Hook, S.E., et al., RNA-Seq analysis of the toxicant-induced transcriptome of the marine diatom, Ceratoneis closterium, Mar. Genomics (2013), http://dx.doi.org/10.1016/j.margen.2013.12.004

275 276 277 278 279 280 281 282 283 284 285 286 287 Q14 288 289

292 293 294 295 296 297 298 299 300 301 302 303 304

4 Table 2 Summary of transcriptome assembly data.

t2:3

Consensus results

t2:4 t2:5 t2:6 t2:7

Number of reads available Number of reads assembled Contig metrics

With rRNA

Without rRNA

t2:8 Q4 t2:9 t2:10

Number of contigs Number of bases Average contig size

12,779 11,579,318 906

4768 5,471,903 1146

320 321 322 323 324 325 Q17 326 327 328 329 330 331 332 333 334 335 336 Q18 337 338

3.4. Differential expression

340 341

When the total number of reads obtained for each cDNA library was compared, considerable variability was seen, with the control reads having approximately 1.7 times the number of reads than that of the copper-exposed library (Table 4). This bias was likely introduced during the cDNA library preparation and has been discussed elsewhere (454 Sequencing Manual, Guidelines for Amplicon Experimental Design, March 2011). Therefore, it is likely that these do not have any ecotoxicological significance. When calculating fold change, we normalized the total number of reads in each library in order to avoid bias by dividing each mapping read by the total number of reads in the treatment. As shown in Fig. 4, the most abundant transcripts in each treatment had a considerable overlap. The most abundant transcripts in the control included a cytochrome p450, an ATP synthase subunit and a cell wall associated partial protein sequence (Supplemental Tables S3–S7). These top three most abundant transcripts were the same across all contaminant treatments. Other abundant transcripts in the cDNA

344 345 346 347 348 349 350 351 352 353 354 355

U

N

339

342 343

F O R O P D

318 319

T

316 Q16 317

C

314 315

E

312 313

R

310 311 Q15

Fig. 1. Distribution of contig lengths for all contigs with rRNA removed (panel A). For reference, distribution of unfiltered read lengths is provided in panel B. Q2

R

308 309

sequences aligned to the diatoms (1242 to the phylum Bacillariophyta and 840 to class Bacillariophyceae). At the family level, 789 contigs aligned best to Phaeodactylaceae, and 380 aligned best to Thalassiosiraceae. To gauge the coverage of the transciptome, all 12,779 contigs (including rRNA) were compared to the set of 458 highly conserved eukaryotic proteins curated by CEGMA (Parra et al., 2008). Approximately 50% (225 of 458) of proteins were identified amongst our contigs using BLASTx (E value b = 1e − 7) (Supplemental Table S2). We would hypothesize that there were comparatively contigs because of differences in transcription levels. It may be that this diatom has a few transcripts that are highly abundant, and the rest are of much lower copy number. Theoretically, the transcripts described in this study could represent the most abundantly expressed only. Additional studies with a much greater sequencing depth may be required to describe the rest of the transcriptome, which may be best carried out on a different sequencing platform. The lack of standardized protocols for data analysis and quality assurance is a major challenge in using next-generation sequencing data for de novo assembly projects (Martin and Wang, 2011; Hornett and Wheat, 2012; Schliesky et al., 2012). In a previous work (Hook et al., in press), we described the importance of quality assurance steps to ensure that contigs are an accurate representation of the transcriptome. Although the current effort identified fewer unique transcripts than was anticipated, the contigs generated in this study were of high quality in that a) they averaged 1146 bp in length when contigs that mapped to rRNA were removed; b) more than 80% of contigs had an ortholog, typically to another diatom, in the BLAST database; c) the portion of the contigs that aligned to the ortholog in the BLAST database was high; and d) a majority of the commonly expressed eukaryotic proteins were included in this sequencing effort. As discussed in a previous work (Hook et al., in press), use of these criteria will increase the probability of obtaining contigs that accurately reflect the sequence composition.

O

306 307

C

305

1,248,932 263,589 (21%)

E

t2:1 t2:2

S.E. Hook et al. / Marine Genomics xxx (2013) xxx–xxx

libraries included a photosystem II protein, pg1 protein, associated protein, and cell wall associated partial protein (Supplemental Tables S3–S7). To ensure that the most sequence information for each treatment was obtained using the 454 pyrosequencing technique, pools of approximately 3 billion individual diatoms cells were used to create cDNA libraries, but no replicate libraries were created for each treatment. Other studies using Roche's 454 pyrosequencing technique have used a similar approach (e.g. Bellin et al., 2009; Hale et al., 2009; Burns et al., 2013). Determining differential expression is also an aim of some of these studies. For instance, other studies have compared expression levels in different tissues of birds and prawns (e.g. Santure et al., 2011; Jung et al., 2012), in different amphipod developmental stages (Zeng et al., 2011), or between stressed and control individual coral colonies (Traylor-Knowles et al., 2011). As it is, the relatively few unique transcripts obtained suggest that there was an insufficient coverage depth to capture rare sequences. The lack of biological replication means that our capacity to determine statistical significance is reduced (Anders and Huber, 2010). Instead, we focused on those contigs with the highest normalized read count. Differentially expressed transcripts were determined by comparing normalized contig abundance in each of the treatments to abundance in control via fold change. Since we only had one cDNA library per

Please cite this article as: Hook, S.E., et al., RNA-Seq analysis of the toxicant-induced transcriptome of the marine diatom, Ceratoneis closterium, Mar. Genomics (2013), http://dx.doi.org/10.1016/j.margen.2013.12.004

356 357 358 359 360 361 362 363 Q19 Q20 364 365 Q21 366 Q22 367 Q23 368 369 Q24 370 371 372 373 Q25 374 375 376 377

S.E. Hook et al. / Marine Genomics xxx (2013) xxx–xxx

5

410 411

3.5.1. Photosynthetic pathways Many of the transcripts identified in this study had roles in photosynthesis. These transcripts included photosystem I proteins (19 contigs), fucoxanthin chlorophyll proteins (27 contigs), ribulose-bisphophate carboxylase oxygenase subunits (66 contigs) and photosystem II proteins (81 contigs) (Supplementary Tables S3–S7). In addition, by mapping to KEGG pathways, 119 contigs were identified as having roles in “carbon fixation in photosynthetic organisms” and 35 contigs had roles in “carbon fixation in prokaryotes” (Supplementary Table 8). Photosynthesis-specific transcripts were abundant in all treatments; in fact, a photosytem II protein was the fourth most abundant contig in the control treatment (Supplementary Table S2). However in comparison to the control, they decreased in abundance after exposure to ammonia, copper and oil (Supplementary Tables S9–S11). Conversely,

430

O

F

pathways. As shown in Table 4, 888 of the original contigs or 812 of the rRNA filtered contigs could be mapped to KEGG pathways. The result of mapping of transcripts to KEGG pathways is shown in Supplementary Table 8. Metabolism related pathways were commonly enriched in all treatments, these included glyoxylate and dicarboxylate metabolism, purine metabolism, thiamine metabolism, pyruvate metabolism, amino and nucleotide sugar metabolism and transamine metabolism. Glycine dehydrogenase (decarboxylating) (EC: 1.4.2), which plays a role in photosynthesis and plant growth was abundant only in simazine and control. The enzyme pathway for porphyrin and chlorophyll metabolism (EC: 1.3.3.3) was abundant only in the simazine treatment. Comparing the functional annotation results to the differential expression results provided the best insight into how physiological pathways were altered in each treatment. These changes are summarized in Table 5 and explored in detail in the next section. However, a substantial number of differentially expressed transcripts were labelled as hypothetical proteins or had no annotation, this included 72 unknown transcripts out of a total 325 transcripts in the ammonia treatment, copper 80 (from total 319), simazine 51 (from total 151), and oil 65 (from total 210).

404 405

409

Of the 12,779 original contigs, 5236 had homologs in the NCBI nr protein database that permitted the assignment of GO annotation. Once rRNA contigs were removed, 3038 of the remaining 4769 contigs could be assigned GO annotation. The distribution of these GO terms is shown in Supplemental Fig. 1. When possible, the contigs were also mapped to KEGG (Kyoto Encyclopaedia of Genes and Genomes)

t3:1 t3:2 t3:3 t3:4 t3:5

Table 3 Summary of the contig annotation (with rRNA contamination removed) by comparison with different databases. Nucleotide orthologs were identified using a BLASTn search in the CSIRO instance of Galaxy, protein orthologs were identified using a BLASTx search in BLAST2GO, both searches had a minimum evalue of 1 × 10−3.

395 396 397 398 399 400

406 407 408

t3:6 t3:7 t3:8 t3:9 t3:10 t3:11

C

393 394

E

391 392

R

389 390

R

387 388

N C O

385 386

U

383 384

P

3.5. Functional annotation

381 382

D

403

379 380

T

401 402

treatment, change was determined conservatively, with those transcripts with fold change greater than 10 and a read count greater than 10 considered differentially expressed. The most abundant contigs in every treatment were highly similar, despite exposure to the chosen contaminants. In addition, the contigs in each treatment with the greatest fold change relative to their expression levels in control were highly similar. Ammonia treated samples had the highest number of differentially expressed transcripts (325), while copper had 319, oil 210 and simazine 151 differentially expressed transcripts. When differences in the numbers of reads for each treatment was accounted for, most of these putatively differentially expressed transcripts were down regulated in all treatments. Comparison of the down regulated genes in the Venn diagrams (Fig. 5) shows a number of differentially expressed transcripts common to all treatments. For instance, there were 132 contigs in common between copper and ammonia exposures 117 in common between oil and ammonia, and 94 in common between copper and oil. This suggests that there may be conserved changes in growth and other physiological properties following exposure to the different contaminants, as discussed in Section 3.5. Alternatively, it may be that this over representation of down regulated genes is a stress response, resulting from the cells favouring homeostasis over growth. This is in contrast to previous microarray based analysis of the transcriptomic response to contaminant exposure in another marine diatom (Osborn and Hook, 2013). There were not a sufficient number of transcripts upregulated to warrant comparison.

E

378

R O

Fig. 2. Proportion of the contig that aligns to a GenBank sequence.

Number of contigs Number of contigs with a nucleotide ortholog Number of contigs with a protein ortholog Number of unique nucleotide orthologs (number of distinct transcripts) Number of contigs with GO Hits Number of contigs with KEGG terms

4769 3131 3997 1660 3038 812

Fig. 3. The taxonomic distribution of the top hit species distribution for the fixed length trimmed MIRA assembly, as generated by metagenomic analysis tools following a BLASTn search. The most common taxonomic groupings are shown in red. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Please cite this article as: Hook, S.E., et al., RNA-Seq analysis of the toxicant-induced transcriptome of the marine diatom, Ceratoneis closterium, Mar. Genomics (2013), http://dx.doi.org/10.1016/j.margen.2013.12.004

412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429

431 432 433 434 435 436 437 438 439 440 441 442 443

453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 471 472

473 474

C

E

R

R

O

C

N

U

t4:1 t4:2 t4:3

3.5.3. Purine metabolism The KEGG pathway to which the greatest number of contigs matched was that of purine metabolism (Supplementary Table S8). In total, 168 contigs with roles in DNA metabolism were identified, including 14 polymerases, 14 helicases and 64 phosphatases. Although there are many potential intracellular purines, one of the most obvious is nucleotide metabolism. Some of these contigs, such as AAA ATPase, were among the most abundant in many treatments (Supplementary treatments 3–7). However, contigs with roles in purine metabolism were also among the most frequently down regulated, particularly following exposure to oil, ammonia, and copper (Supplementary Tables S9–S11). Previous studies with this alga have shown that contaminant exposure can inhibit cell division while having no effect on photosynthesis, respiration, ATP production and the electron transport systems (Stauber and Forence, 1987). Other prior studies have also found similar changes in gene expression following contaminant exposure. DNA replication and cell cycle related genes have been reported to be differentially expressed after copper toxicity in the green alga C. reinhardtii (Jamers et al., 2006) and similar to our study the purine metabolism KEGG pathway was abundant in rice (Oryza sativa) gene expression studies exposed to atrazine (Zhang et al., 2012).

500 501

3.5.4. Nutrient uptake and metabolism A number of transcripts with roles in nutrient uptake and metabolism were also identified, including two nitrate transporters and four others identified via GO mapping with roles in nitrate transport. Two silicon transporters were also identified. Silica deposition is important for diatom cell walls. Sil3, which is involved in this process, was repressed in the marine diatom T. pseudonana upon exposure to PAHs, a major component of oil (Bopp and Lettieri, 2007; Carvalho et al., 2011). One of the nitrate transporters was less common in the Cu, oil, and ammonia exposed cDNA libraries (Supplemental Tables S9–S11).

521

Table 4 Summary of mapping statistics for individual treatments, once ribosomal RNA contamination is removed. The number of reads, after filtering and trimming, is provided for each treatment. The number of those filtered reads that could be aligned to the assembled contigs if also provided.

t4:4 t4:5 t4:6

Fig. 4. Overlap amongst the most abundant contigs in each treatment. Since there are five treatments, two figures are presented: the upper compares control, Cu, oil and simazine, the lower panel compares control, Cu, ammonia and simazine.

T

3.5.2. Respiration A large number of the contigs identified in this study had roles in 475 the metabolic pathways of respiration. As shown in Supplementary 476 Table S8, 103 contigs with roles in the KEGG pathway “Glyoxylate and 477 dicarboxylate metabolism” were identified, along with 72 contigs with 478 roles in “pyruvate metabolism”, 66 contigs with roles in “oxidative 479 phosphorylation” and 54 with roles in “gycolysis/gluconeogenesis.” 480 Despite their physiological importance, these respiration related 481 contigs were not among the most abundant in any treatment (Supple482 mentary Tables S3–S7). However, these contigs were among the most 483 frequently decreased in abundance following contaminant exposure. 484 For instance, contigs encoding fructose-bisphosphate aldolase, which 485 have roles both in carbon fixation and in glycolysis, were decreased 486 in abundance by as much as ten fold following exposure to ammonia, 487 copper or oil (Supplementary Tables S9–S11). Other sequences involved 488 in glycolysis, including 6-phosphofructokinase 3, pyruvate kinase and 489 pyruvate dehydrogenase, were at least ten fold less abundant in the 490 ammonia, oil or simazine exposed cDNA (Supplementary Tables S9, 491 Q26 S11, S12). Other contigs with putative roles in pyruvate metabolism 492 such as formate acetyltransferase, pyruvate dehydrogenase and 493 pyruvate carboxylase were also at least ten times less abundant in the 494 contaminant exposed cDNA libraries (Supplementary Tables S9–S12). 495 Pyruvate- and glycolysis-related genes have also been observed to be 496 down regulated after copper toxicity in the algae Chlamydomonas 497 Q27 reinhardtii (Jamers et al., 2006) and Scytosiphon gracilis (Contreras 498 et al., 2010). Taken together, these data suggest that respiration is 499 substantially altered by contaminant exposure (Table 5).

F

451 452

O

450

R O

448 449

P

446 447

exposure to simazine resulted in increased abundance of photosynthetic related transcripts compared to the control (Table 5; Supplementary Table S12). A reduction in photosynthesis could explain the decrease in cell growth observed when algae were exposed to sub-lethal concentrations for these contaminants. Previous studies observed a decrease in photosynthetic rates in response to stressors and suggested it as a specific action to protect against oxidative stress (Fortes et al., 2008; Qian et al., 2009; Saibo et al., 2009). Specifically, exposure to oil or PAHs has commonly been reported to decrease photosynthesis activity in microalgae (Carrera-Martinez et al., 2010; Romero-Lopez et al., 2012) and down regulate photosynthesis related genes such as light harvesting proteins (Bopp and Lettieri, 2007). Interestingly, a small increase in photosynthetic transcript abundance was unique to the simazine treatment (Supplementary Table S12). This may be a specific response to the mechanism of toxicity for this herbicide since it targets photosystem II (Magnusson et al., 2010). Simazine exposure resulted in the presence of more chloroplast-localised transcripts, specifically, cell wall associated transcripts as well as coproporphyrinogen III oxidase, a transcript involved in the biosynthetic pathway of chlorophyll (Ishikawa et al., 2001). Its corresponding enzyme pathway, porphyrin and chlorophyll metabolism was also only observed in this treatment (Supplementary Table S12). Glycine dehydrogenase was abundant in both control and simazine treatments and is involved in controlling photosynthesis and plant growth (Timm et al., 2012). Studies of triazinebased herbicide exposure on plant gene expression have demonstrated both increased and decreased abundance of photosynthetic genes, especially those associated with light reactions (Qian et al., 2008; Zhu et al., 2009).

D

444 445

S.E. Hook et al. / Marine Genomics xxx (2013) xxx–xxx

E

6

Number of reads Total number of alignments

Control

Ammonia

Cu

Oil

Simazine

317,218 295,176 (93%)

260,918 255,722 (98%)

186,799 183,184 (98%)

253,022 245,858 (97%)

204,905 166,438 (81%)

Please cite this article as: Hook, S.E., et al., RNA-Seq analysis of the toxicant-induced transcriptome of the marine diatom, Ceratoneis closterium, Mar. Genomics (2013), http://dx.doi.org/10.1016/j.margen.2013.12.004

502 503 504 505 506 507 508 509 510 511 512 513 514 Q28 515 516 517 518 Q29 519 520

522 523 524 525 526 527 528 529 530

S.E. Hook et al. / Marine Genomics xxx (2013) xxx–xxx

7

540 541

4. Conclusions

577

This paper describes a partial transcriptome of the benthic diatom Ceratoneis closterium and provides new insight into a taxonomic group that has been underrepresented in the genomics literature, despite significant ecological importance. This sequencing effort describes 4768 non-rRNA contigs. Although the 1660 unique transcripts described via our sequencing project are not the full complement of expressed genes, this effort will provide valuable information regarding not only toxicant and other stressor responses, but other important physiological process including photosynthesis and respiration. The data also enhances the ecotoxicogenomic resources available for this ecologically important taxonomic group.

578 579

5. Uncited references

589 Q32

Fig. 5. Overlap amongst the most differentially expressed contigs in each treatment. Panel A shows the degree of overlap among the most 50 upregulated transcripts for each treatment; panel B depicts the overlap amongst the 50 most down regulated contigs in each treatment.

N C O

Q3

R

R

E

C

T

E

D

P

R O

O

F

p450 which have a known role in the removal of xenobiotic substances (Gonzalez, 2005). Senescence-associated proteins were also abundant and are involved in programmed cell death (Pennell and Lamb, 1997). These senescence-associated transcripts have been reported in previous stress related studies including exposure to abiotic stressors in chickpeas (Mantri et al., 2007), nitrogen starvation in the green algae, Micractinium pusillum (Li et al., 2012) and salinity stress in the crop cultivar Setaria italica (Puranika et al., 2011). In the simazine treatment only, senescence-associated proteins had decreased abundance and transcripts with hydrolase activity had increased abundance relative to all other treatments and control (Supplementary Table S12). Another stress indicator, serine palmitoyltransferase, which plays a role in apoptosis (Hanada, 2004), was down regulated in ammonia and copper treatments. Other potential stress related transcripts included serine carboxypeptidase which was present in the copper and simazine treatments, as well as the control. The variety of roles this protease enzyme may have is extensive, though it has been reported to increase protein turnover during wound stress and also increase in activity during senescence (Jiang et al., 1999; Schaller, 2004). The high abundance of mitochondrial alternative oxidase transcripts in the copper and ammonia treatments is also linked to toxicity as it is thought to be influenced by stress stimuli (Vanlerberghe and McIntosh, 1997). In addition, the antioxidant glutathione peroxidase (Mallick and Mohn, 2000) was down regulated in response to both ammonia and oil. We were unable to discern trends related to the modes of action of the different contaminants. The lack of treatment-specific transcripts may be due to the algae using stress responses not identified by this technique. Alternatively it may be that our study used contaminant concentrations that were too high — resulting in generalised stress responses, which, as stated previously, are in contrast to our previous studies of toxicant related changes in the transcriptome of marine diatoms (Osborn and Hook, 2013). Unexpectedly, the majority of transcripts differentially expressed in all the treatments were down regulated, perhaps as a result of the slightly reduced growth rate in these cells. Compounding the difficulty in identifying trends in this study was the many unknown or hypothetical proteins.

531

Nitrate transporters have also been observed in P. tricornutum exposed to ammonia (Osborn and Hook, 2013) and decreased expression of 533 genes controlling nitrogen assimilation have been observed upon expo534 sure to heavy metals and atrazine in the cyanobacteria Microcystis 535 Q30 aeruginosa (Qian et al., 2012). 536

U

532

539

3.5.5. Stress response genes There were a number of transcripts with involvement in stress responses, an expected response to compensate for contaminant exposure, including ABC transporters, heat shock proteins and cytochrome

t5:1 t5:2 t5:3

Table 5 Summary of the changes in transcript abundance measured in each treatment, as compared to controls. Abundant transcripts are among the 50 most frequently counted contigs; whereas differentially expressed contigs have a read count greater than 10 and fold change greater than 10.

537 538

Akerfelt et al., 2010 Mattick and Makunin, 2006

t5:4

Functional group

Ammonia

Copper

Oil

Simazine

t5:5 t5:6 t5:7 t5:8 t5:9

Photosynthesis Respiration Purine metabolism Nutrient cycling Stress response

Decreased abundance Down regulated Down regulated Down regulated

Decreased abundance Down regulated Down regulated Down regulated

Decreased abundance Down regulated Down regulated Down regulated

Increased abundance Down regulated Abundant Increased abundance (senescence related proteins)

Please cite this article as: Hook, S.E., et al., RNA-Seq analysis of the toxicant-induced transcriptome of the marine diatom, Ceratoneis closterium, Mar. Genomics (2013), http://dx.doi.org/10.1016/j.margen.2013.12.004

542 543 544 545 546 547 548 549 550 551 552 Q31 553 554 555 556 557 558 559 560 561 562 563 564 565 566 567 568 569 570 571 572 573 574 575 576

580 581 582 583 584 585 586 587 588

590 591

610

Appendix A. Supplementary data

611 612

Supplementary data to this article can be found online at http://dx. doi.org/10.1016/j.margen.2013.12.004.

613

References

614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654 655 656 657 658 659 660 661 662

Adams, M.S., Stauber, J.L., Binet, M.T., Molloy, R., Gregory, D., 2008. Toxicity of a secondary-treated sewage effluent to marine biota in Bass Strait, Australia: development of action trigger values for a toxicity monitoring program. Mar. Pollut. Bull. 57, 587–598. Akerfelt, M., Morimoto, R., Sistonen, L., 2010. Heat shock factors: integrators of cell stress, development and lifespan. Nat. Rev. Mol. Cell Biol. 11, 545–555. Aurand, D., Coelho, G. (Eds.), 2005. Cooperative Aquatic Toxicity Testing of Dispersed Oil and the Chemical Response to Oil Spills: Ecological Effects Research Forum (CROSERF). Technical Report 07-03. Ecosystem Management and Associates, Inc., Lusby, MD (105 pp. + Appendices). Binet, M.T., Adams, M.S., Stauber, J.L., King, C.K., Doyle, C.J., Lim, R.P., Laginestra, E., 2001. Toxicity assessment of leachates from Homebush Bay landfills. Australas. J. Ecotoxicol. 9, 7–18. Blankenberg, D., Gordon, A., Von Kuster, G., Coraor, N., Taylor, J., Nekrutenko, A., The Galaxy Team, 2010. Manipulation of FASTQ data with Galaxy. Bioinformatics 26, 1783–1785. Bopp, S.K., Lettieri, T., 2007. Gene regulation in the marine diatom Thalassiosira pseudonana upon exposure to polycyclic aromatic hydrocarbons (PAHs). Gene 396, 293–302. Burns, G., Thorndyke, M.C., Peck, L.S., Clark, M.S., 2013. Transcriptome pyrosequencing of the Antarctic brittle star Ophionotus victoriae. Mar. Genomics 9, 9–15. Carrera-Martinez, D., Mateos-Sanz, A., Lopez-Rodas, V., Costas, E., 2010. Microalgae response to petroleum spill: an experimental model analysing physiological and genetic response of Dunaliella tertiolecta (Chlorophyceae) to oil samples from the tanker Prestige. Aquat. Toxicol. 97, 151–159. Carvalho, R.N., Burchardt, A.D., Sena, F., Mariani, G., Mueller, A., Bopp, S.K., Umlauf, G., Lettieri, T., 2011. Gene biomarkers in diatom Thalassiosira pseudonana exposed to polycyclic aromatic hydrocarbons from contaminated marine surface sediments. Aquat. Toxicol. 101, 244–253. Chevreux, B., Wetter, T., Suhai, S., 1999. Genome sequence assembly using trace signals and additional sequence information. Computer Science and Biology: Proceedings of the German Conference on Bioinformatics (GCB), 99, pp. 45–56. Conesa, A., Götz, 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. Contreras, L., Moenne, A., Gaillard, F., Potin, P., Correa, J.A., 2010. Proteomic analysis and identification of copper stress regulated proteins in the marine alga Scytosiphon gracilis (Phaeophyceae). Aquat. Toxicol. 96, 85–89. Craft, J.A., Gilbert, J.A., Temperton, B., Dempsey, K.E., Ashelford, K., Tiwari, B., Hutchinson, T.H., Chipman, J.K., 2010. Pyrosequencing of Mytilius galloprovincialis cDNAs: tissuespecific expression patterns. PLoS One 5, e8875. Falkowski, P.G., Katz, M.E., Knoll, A.H., Quigg, A., Raven, J.A., Schofield, O., Taylor, F.J.R., 2004. The evolution of modern eukaryotic phytoplankton. Science 305, 354–360. Fortes, et al., 2008. Organogenic nodule development in hop (Humulus lupulus L.): transcript and metabolic responses. BMC Genomics 9, 445. Franklin, N.M., Stauber, J.L., Adams, M.S., 2005. In: Ostrander, G.K. (Ed.), Techniques in Aquatic Toxicology, volume 2. CRC Press, Boca, Ranta, FL, USA, pp. 735–756. Giardine, B., Riemer, C., Hardison, R.C., Burhans, R., Elnitski, L., Shah, P., Zhang, Y., Blankenberg, D., Albert, I., Taylor, J., Miller, W., Kent, W.J., Nekrutenko, A., 2005.

C

E

R

608

R

606 607

O

604 605

C

602 603

N

600 601

U

598 599

F

609

The Systems Biology Initiative acknowledges support from the EIF Super Science Scheme, the NSW State Government Science Leveraging Fund and the University of New South Wales. Additional support was made available by the CSIRO Wealth from Oceans Flagship and the CSIRO Bioinformatics Core. The authors acknowledge the assistance of Ms Sarah Stevenson (CSIRO, Land and Water) with cell manipulation and exposures. Mr. Jason Koval (Ramaciotti Centre for Gene Function Analysis, University of New South Wales) performed the 454 sequencing. The primary data and workflows can be found at: Sharon Hook; Hannah Osborn; Francesca Gissi; Philippe Moncuquet; Natalie Twine; Marc Wilkins; et al. (2013): Raw reads and final assembled contigs for pennate diatom project. v1. CSIRO. Data collection http://dx.doi.org/ 10.4225/08/525BCD6E6F145. This manuscript was improved by CSIRO internal reviewers, Sarah Stevenson and Graeme Batley, and by two anonymous reviewers

596 597

O

595

R O

Acknowledgements

T

594

Galaxy: a platform for interactive large-scale genome analysis. Genome Res. 15, 1451–1455. Goecks, J., Nekrutenko, A., Taylor, J., The Galaxy Team, 2010. Galaxy: a comprehensive approach for supporting accessible, reproducible, and transparent computational research in the life sciences. Genome Biol. 11, R86 (25). Gonzalez, F., 2005. Role of cytochromes P450 in chemical toxicity and oxidative stress: studies with CYP2E1. Mutat. Res. Fundam. Mol. Mech. Mutagen. 569 (1–2), 101–110. Guillard, R.R.L., Ryther, J.H., 1962. Studies on marine planktonic diatoms. I. Cyclotella nana Hustedt and Detonula confervacea (Cleve) Gran. Can. J. Microbiol. 8, 229–239. Hanada, K., 2004. Serine palmitoyltransferase, a key enzyme of sphingolipid metabolism. Biochim. Biophys. Acta 1632 (1–3), 16–30. Hogan, A.C., Stauber, J.L., Pablo, F., Adams, M.S., Lim, R.P., 2005. The development of marine toxicity identification evaluation (TIE) procedures using a unicellular alga Nitzschia closterium. Arch. Environ. Contam. Toxicol. 48, 433–443. Hook, S.E., Osborn, H.L., 2012. Comparison of toxicity and transcriptomic profiles in a diatom exposed to oil, dispersants, dispersed oil. Aquat. Toxicol. 124/125, 139–151. Hook, S.E., Lampi, M.A., Febbo, E.J., Ward, J.A., Parkerton, T.F., 2010. Temporal patterns in the transcriptomic response of rainbow trout, Oncorhynchus mykiss, to crude oil. Aquat. Toxicol. 99, 320–329. Hook, S.E., Twine, N.A., Simpson, S.L., Spadaro, D.A., Moncuquet, P., Wilkins, M., 2013. Toxicant-induced changes in the transcriptome of the amphipod Melita plumulosa. Aquat. Toxicol. (in press). Hornett, E.A., Wheat, C.W., 2012. Quantitative RNA-Seq analysis in non-model species: assessing transcriptome assemblies as a scaffold and the utility of evolutionary divergent genomic reference species. BMC Genomics 13, 361. Hudson, M.E., 2008. Sequencing breakthroughs for genomic ecology and evolutionary biology. Mol. Ecol. Resour. 8, 3–17. Ishikawa, A., Okamoto, H., Iwasaki, Y., Asahi, T., 2001. A deficiency of coproporphyrinogen III oxidase casues lesion formation in Arabidopsis. Plant J. 27 (2), 89–99. Jiang, W.B., Lers, A., Lomaniec, E., Aharoni, N., 1999. Senescence-related serine protease in parsley. Phytochemistry 50, 377–382. Levy, J.L., Stauber, J.L., Jolley, D.F., 2007. Sensitivity of marine microalgae to copper: the effect of biotic factors on copper adsorption and toxicity. Sci. Total Environ. 387, 141–154. Li, Y., Fei, X., Deng, X., 2012. Novel molecular insights into nitrogen starvation-induced triaclyglycerols accumulation revealed by differential gene expression analysis in green algae Micractinium pusillum. Biomass Bioenergy 42, 199–211. Lowe, C.D., Mello, L.V., Samatar, N., Martin, L.E., Montagnes, D.J.S., Watts, P.C., 2011. The transcriptome of the novel dinoflagellate Oxyrrhis marina (Alveolata: Dinophycae): response to salinity examined by 454 sequencing. BMC Genomics 12, 519. Magnusson, M., Heimann, K., Quayle, P., Negri, A.P., 2010. Additive toxicity of herbicide mixtures and comparative sensitivity of tropical benthic microalgae. Mar. Pollut. Bull. 60, 1978–1987. Mallick, N., Mohn, F., 2000. Reactive oxygen species: response of algal cells. J. Plant Physiol. 157, 183–193. Mantri, N.L., Ford, R., Coram, T.E., Pang, C.K., 2007. Transcriptional profiling of chickpea genes differentially regulated in response to high salinity, cold and drought. BMC Genomics 8, 303. Martin, J.A., Wang, Z., 2011. Next-generation transcriptome assembly. Nat. Rev. Genet. 12, 671–682. Mattick, J.S., Makunin, I.V., 2006. Non-coding RNA. Hum. Mol. Genet. 15, R17–R29. McCarthy, F.M., Wang, N., Magee, G.B., Nanduri, B., Lawrence, M.L., Camon, E.B., Barrell, D.G., Hill, D.P., Dolan, M.E., Williams, W.P., Luthe, D.S., Bridges, S.M., Burgess, S.C., 2006. AgBase: a functional genomics resource for agriculture. BMC Genomics 7, 229. Moreno-Garrido, I., Hampel, M., Lubián, L.M., Blasco, J., 2003. Sediment toxicity tests using benthic marine microalgae Cylindrotheca closterium (Ehrenberg) Lewin and Reimann (Bacillariophyceae). Ecotoxicol. Environ. Saf. 54, 290–295. Osborn, H.L., Hook, S.E., 2013. Using transcriptomic profiles in the diatom Phaeodactylum tricornutum to identify and prioritize stressors. Aquat. Toxicol. 138–139, 12–25. Parker, M.S., Mock, T., Armbrust, E.V., 2008. Genomic insights into marine microalgae. Annu. Rev. Genet. 42, 619–645. Pennell, R., Lamb, C., 1997. Programmed cell death in plants. Plant Cell 9, 1157–1168. Plaxton, W., Podesta, F., 2006. The functional organization and control of plant respiration. Crit. Rev. Plant Sci. 25, 159–198. Poynton, H.C., Vulpe, C.D., 2009. Ecotoxicogenomics: emerging technologies for emerging contaminants. J. Am. Water Resour. Assoc. 45, 83–96. 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. Puranika, S., Jha, S., Srivastava, P.S., Sreenivasulu, N., Prasad, M., 2011. Comparative transcriptome analysis of contrasting foxtail millet cultivars in response to shortterm salinity stress. J. Plant Physiol. 168, 280–287. Qian, H.F., Sheng, G.D., Liu, W.P., Lu, Y.C., Liu, Z.G., Fu, Z.W., 2008. Inhibitory effects of atrazine on Chlorella vulgaris as assessed by real-time polymerase chain reaction. Environ. Toxicol. Chem. 27, 182–187. Qian, H., Li, J., Sun, L., Chen, W., Sheng, G., Liu, W., Fu, Z., 2009. Combined effect of copper and cadmium on Chlorella vulgaris growth and photosynthesis-related gene transcription. Aquat. Toxicol. 94, 56–61. Qian, H., Pan, X., Chen, J., Zhou, D., Chen, Z., Zhang, L., Fu, Z., 2012. Analyses of gene expression and physiological changes in Microcystis aeruginosa reveal the phytotoxicities of three environmental pollutants. Ecotoxicology 21, 847–859. Quast, C., Pruesse, E., Yilmaz, P., Gerken, J., Schweer, T., Yarza, P., Peplies, J., Glöckner, F.O., 2013. The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. Nucleic Acids Res. 41 (D1), D590–D596.

P

Plaxton and Podesta, 2006 Thompson et al., 1988

D

592 593

S.E. Hook et al. / Marine Genomics xxx (2013) xxx–xxx

E

8

Please cite this article as: Hook, S.E., et al., RNA-Seq analysis of the toxicant-induced transcriptome of the marine diatom, Ceratoneis closterium, Mar. Genomics (2013), http://dx.doi.org/10.1016/j.margen.2013.12.004

663 664 665 666 667 668 669 670 671 672 673 674 675 676 677 678 679 680 681 682 683 684 Q34 685 686 687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 729 730 731 732 733 734 735 736 737 738 739 740 741 742 743 744 745 746 747 748

S.E. Hook et al. / Marine Genomics xxx (2013) xxx–xxx Romero-Lopez, J., Lopez-Rodas, V., Costas, E., 2012. Estimating the capability of microalgae to physiological acclimatization and genetic adaptation to petroleum and diesel oil contamination. Aquat. Toxicol. 124–125, 227–237. Saibo, N.J., Lourenco, T., Oliveira, M., 2009. Transcription factors and regulation of photosynthetic and related metabolism under environmental stresses. Ann. Bot. 103, 609–623. Schaller, A., 2004. A cut above the rest: the regulatory function of plant proteases. Planta 220, 183–197. Schliesky, S., Gowik, U., Weger, A.P.M., Bautigam, A., 2012. RNA-Seq assembly — are we there yet? Front. Plant Sci. 3, 220. Stauber, J.L., Forence, T.M., 1987. Mechanism of toxicity of ionic copper and copper complexes to algae. Mar. Biol. 94, 511–519. Stauber, J.L., Tsai, J., Vaughan, G., Peterson, S.M., Brockbank, C.I., 1994. Algae as indicators of toxicity of BEKM effluents. National Pulp Mills Research Program Technical Report Series No. 3. CSIRO, Canberra, Australia, pp. 1–83. Thompson, A.S., Rhodes, J.C., Pettman, I., 1988. Culture Collection of Algae and Protozoa: Catalogue of Strains. Natural Environmental Research Council, Swindon, UK.

Timm, S., Florian, A., Arrivault, S., Stitt, M., Fernie, A., Bauwe, H., 2012. Glycine decarboxylase controls photosynthesis and plant growth. FEBS Lett. 586, 3692–3697. Vanlerberghe, G.C., McIntosh, L., 1997. Alternative oxidase: from gene to function. Annu. Rev. Plant Physiol. Plant Mol. Biol. 48, 703–734. Villeneuve, D.L., Garcia-Reyero, N., 2011. Vision & strategy predictive ecotoxicology in the 21st century. Environ. Toxicol. Chem. 30, 1–8. Zeng, V., Villanueva, K.E., Ewen-Campen, B.S., Alwes, F., Browne, W.E., Extavour, C.G., 2011. De novo assembly and characterization of a maternal and developmental transcriptome for the emerging model crustacean Parhyale hawaiensis. BMC Genomics 12, 581. Zhang, J.J., Zhou, Z.S., Song, J.B., Liu, Z.P., Yang, H., 2012. Molecular dissection of atrazineresponsive transcriptome and gene networks in rice by high-throughput sequencing. J. Hazard. Mater. 219–220, 57–68. Zhu, J., Patzoldt, W.L., Radwan, O., Tranel, P.J., Clough, S.J., 2009. Effects of photosystem-ii-interfering herbicides atrazine and bentazon on the soybean transcriptome. Plant Genome 2, 191–205.

F

749 750 751 752 753 754 755 756 757 758 759 760 761 762 763 764 765

9

U

N C O

R

R

E

C

T

E

D

P

R O

O

784

Please cite this article as: Hook, S.E., et al., RNA-Seq analysis of the toxicant-induced transcriptome of the marine diatom, Ceratoneis closterium, Mar. Genomics (2013), http://dx.doi.org/10.1016/j.margen.2013.12.004

766 767 768 769 770 771 772 773 774 775 776 777 778 779 780 781 782 783