BACKGROUND
the ecologically relevant soil-living worm enchytraeus crypticus  has been widely used as model species in soil ecotoxicology for over two decades, both under international organization for standardization  and organization for economic co-operation and development  guidelines  <cit> .
ecotoxicological tests generally focus on survival and reproduction endpoints, i.e.
organism effects, whereas in ecotoxicogenomics, responses at the transcriptome level are studied aiming to reconstitute the mechanisms preceding effects at organismal and/or population level .
in-depth information about the molecular mechanisms is provided allowing a better understanding of molecular pathways activated by toxicity, which will not be possible when using the classical tests alone  <cit> .
previously, transcript information with classical sanger sequencing was generated for the test species enchytraeus albidus <cit>  and has been very useful to further explore the response mechanisms of enchytraeids to various stressors, such as pesticides  <cit> , heavy metals, nanoparticles  <cit> , and soil properties  <cit> .
however, the drawback of the current e. albidus microarray is that it contains a low number of transcripts, so a full transcriptomic assessment of stress response cannot be executed for this model.
thanks to the progress in sequencing technology, transcriptome and genome information can now be obtained for any non-genomic model species  <cit> .
transcriptome sequencing is of special interest, since it is directly targeted to obtain sequence information from transcribed regions in the genome  <cit> .
to pursue a comprehensive transcriptome of an enchytraeid test species we selected e. crypticus instead of the more widely used e. albidus, despite the seemingly lower sensitivity to toxicants  <cit>  of this species.
several advantages over the use of e. albidus have been recently identified, such as easy culturing conditions in agar media, higher reproductive rate, lack of cryptic speciation, and a shorter reproduction cycle   <cit> .
given the clear advantages of using e. crypticus for soil quality assessment, the study presented here aims at studying its transcriptome by de novo transcriptome pyrosequencing and annotation.
a previous transcriptome study showed that it is beneficial to pool rna from multiple treatments in order to obtain a wider range of expressed transcript sequences  <cit> .
thus, we decided to set up a similar approach to get transcriptome-wide information from e. crypticus and sampled mrna from multiple exposure conditions, including chemical  <cit> , physicochemical  and developmental stages.
the quality of the assembled transcriptome was critically assessed according to the guideline provided by gibson et al.
<cit>  by evaluating sequence quality, assembly quality, overall transcriptome completeness, and transcriptional activity.
subsequently, a high-density microarray platform was constructed and validated using zinc-exposed organisms.
this platform will be a valuable resource to conduct gene expression studies in future ecotoxicogenomic research.
methods
test organism culturing conditions and experiment overview
the annelid enchytraeus crypticus , commonly known as potworm, has been cultured for several years in laboratory at the vu university amsterdam.
cultures are kept in agar substrate prepared with aqueous soil extract, fed ad libitum with ground oatmeal, at 16°c with 75% relative humidity and 16/8 h light/dark photoperiod  <cit> .
for this study we selected, for each biological replicate per exposure, 5– <dig> enchytraeus crypticus adults with developed clitella, similar size, and an average fresh weight of  <dig>   mg per adult.
in total,  <dig> different experimental conditions were included:  <dig> chemical exposures  <dig> temperature treatments,  <dig> soil moisture treatments,  <dig> ph treatments,  <dig> different developmental stages, plus associated control conditions .
test soils
the reference natural lufa  <dig>  soil  was used in all tests, except for the ph treatments, which were performed using oecd artificial soil  <cit> .
lufa  <dig>  soil is a loamy sand soil sieved at 2 mm with ph   <dig> , cation exchange capacity  10 meq/100 g, maximum water-holding capacity  of 46%, organic matter content of 4%, grain size distribution of 6% clay, 14% silt and 80% sand.
for testing, control moisture was adjusted to 50% of the whcmax.
the oecd artificial soil was prepared according to the standard procedure  <cit>  and further used to test various phs adjusting it by adding predetermined amounts of calcium carbonate .
treatments
toxicity tests were performed following the standard guidelines  <cit> .
chemical test concentrations were based on previous results  <cit>  from which the 3-week reproduction ec <dig> was estimated.
this ec <dig> was used as a reference exposure concentration in the short-term  experiments from which the mrna pool was derived.
test chemicals  included cadmium chloride hemipentahydrate , carbendazim , phenanthrene , pentachloroaniline  and  <dig> -dichloroaniline .
chemical formulas, chemical abstract service  numbers and concentrations are described in additional file  <dig>  non-water soluble chemicals were solubilized in acetone.
for the acetone control and acetone chemical spiking the following was performed:  <dig> gr of dry lufa  <dig>  soil was spiked with acetone  or the appropriate organic-acetone solution for phe, car, pca and dca, equilibrated for 24 h, left to evaporate overnight, mixed with dry soil  and homogeneously moistened to 50% of the whcmax.
the soil batches were used immediately after being prepared.
for the metal exposures , dry lufa  <dig>  soil  was spiked with the specific metal-aqueous solution  to reach the nominal concentration and 50% of the whcmax in the soil batch.
the soil was kept under standard conditions while equilibrating for 21 days prior to exposure.
two replicates were prepared per test condition .
per replicate, twenty adult e. crypticus were introduced into each glass vial  containing 30 g of soil prepared as described above.
food  was supplied and vessels were closed with perforated aluminum foil.
the exposures lasted 48 h at 20°c, with 75% relative humidity and 16/8 h light/dark photoperiod.
after 48 h of exposure, the individuals in each replicate were sampled, carefully rinsed in deionised water, transferred to a cryotube, snap-frozen with liquid nitrogen and stored at -80°c.
zinc  exposures were performed in order to validate the microarray platform.
the ec <dig> for inhibition of reproduction was used .
for the microarray analysis, e. crypticus adults were exposed to 145 mg zn/kg soil  and control soil for 4 days using four replicates, containing 5– <dig> individuals per replicate.
temperature, moisture and ph treatments were performed according to the guidelines, with standard test conditions as described in the previous paragraph, except for the specific factor being studied .
the temperature exposures lasted  <dig> ,  <dig>   <dig> or 96 h, at the temperatures - <dig>   <dig>   <dig>   <dig>   <dig> and 40°c.
the moisture exposures lasted  <dig>   <dig> or 96 h, using soils moistened to  <dig>   <dig> and 90% of whcmax.
the ph exposures lasted 48 h, using soils with ph  <dig>   <dig>  and  <dig>
enchytraeus crypticus samples in distinct developmental stages were taken from synchronization of cultures, which were obtained by transferring thirty adults with developed clitella to fresh culture boxes that were kept at identical culturing conditions as the source culture.
every 6 hours, cocoons were collected , immediately snap-frozen with liquid nitrogen and stored at -80°c.
twenty synchronized cocoons were carefully transferred into fresh culture media and kept at 20°c.
after twenty days twenty juveniles were collected into a cryotube, immediately snap-frozen with liquid nitrogen and stored at -80°c.
rna isolation, rna pool and cdna library construction
total rna was collected from the  <dig> individual samples of enchytraeus crypticus as listed in additional file  <dig>  all rna extractions were performed using the sv total rna isolation system  following manufacturer’s protocol.
all samples were homogenized with liquid nitrogen in the cryotube before rna extraction.
addition of sterilized quartz sand to the liquid nitrogen was needed to disrupt and homogenize cocoons sufficiently.
the overall quality of rna samples was evaluated in terms of purity and integrity of rna by means of a nanodrop nd- <dig> uv–vis spectrometer  and agarose gel electrophoresis.
rna sample quality was verified regarding high rna concentration, absorbance ratios a260/a <dig> in the range  <dig>  -  <dig> , and a260/a <dig> above  <dig> .
samples with lower absorbance ratio were ethanol-precipitated in order to improve the quality.
equivalent amounts of rna mass per test condition were pooled together, with a total of 10 μg rna from all samples of e. crypticus.
figure  <dig> illustrates the proportion of rna input per test condition represented in the pool that was used for sequencing.
normalized 454-sequencing libraries were constructed from an equal-molar pool of rna obtained from the unique exposure samples described above using the procedures described in meyer et al.
<cit> .
after the final purification step the library was stored at -20°c until sequencing.
transcriptome sequencing and assembly
cytochrome oxidase i and 18s rrna gene sequencing of  <dig> individuals yielded only one polymorphic site , providing strong evidence that this strain is highly inbred.
thus, the assembly procedure is not hampered by the presence of genetic variation within the cdna library.
this library was submitted to one full-plate sequencing run using a  <dig> roche gs flx pyrosequencing instrument using titanium chemistry , following manufacturer’s protocol and methods previously described  <cit> .
after  <dig> sequencing, the generated sequence reads were assembled using newbler  in de novo mode and default parameters.
transcriptome annotation
transcriptome annotation was performed through the isga transcriptome analysis pipeline, available from the center for genomics and bioinformatics at indiana university .
first, sequence homology to known metazoan proteins was obtained by submitting contigs and singletons to blastx searches against ncbi’s non-redundant database and dbest  <cit> .
moreover, protein domains were identified among the  <dig> frame translations of the assembled sequences using pfam, tigrfam and hmmer <dig> searches.
open reading frames were determined with the orfpredictor software on the proteomics server of youngstown state university  <cit> .
orthology and paralogy was assigned by blastx against orthomcl databases  <cit> .
assembly quality was tested by retrieving blastx hits against the drosophila orthologs in the cegma core eukaryotic genes dataset  <cit> .
transcriptome quality assessment was performed according to the guidelines described by gibson et al.
<cit> .
in short, open reading frame  lengths, est nucleotide lengths and gc contents of annotated and unannoted ests were compared using a welch’s t-test in r version  <dig> .
furthermore, transcript abundance  was compared between contigs and singletons using welch’s t-test.
finally, the interactive pathway explorer ipath was used to identify metabolic, regulatory and biosynthesis pathways based on enzyme code  numbers  <cit> .
probe design for microarray construction
oligonucleotide probes were designed according to the thermodynamics/technical features of gene expression microarrays within the earray webtool of agilent .
the workflow is summarized into these consecutive steps: the target sequences were selected by length threshold ≥ <dig> bases, converted into sense orientation, and then submitted to earray for probe design.
in order to avoid cross-hybridization, the probes that could cross-hybridize with any sequence other than the specific target , were discarded.
a final microarray design was created from the probe group  that were printed in duplicate  onto a  <dig> x <dig> k platform.
the array design was supplemented with  <dig>  additional agilent control probes , so that in total  <dig>  probes were present.
microarray validation and data analysis
rna samples were extracted from control and ec <dig> zn-exposed e. crypticus adults.
the rna was then precipitated with ethanol to remove any guanidine thiocyanate traces from the sv rna extraction procedure and to increase its concentration.
rna concentration and quality was determined on the nanodrop nd <dig> spectrophotometer  and bioanalyzer  <dig> .
rna samples were labeled using the low input quick amp kit  according to the manufacturer’s instructions and as described  <cit> .
four replicates were labeled using cy-3/cy- <dig> dye swaps.
labeled samples were hybridized to the microarray as described above.
microarrays were scanned using an agilent microarray scanner and raw data was extracted with feature extraction software version  <dig> .
<dig> .
microarray raw intensities were loess normalized with a background subtraction according to the edwards method which contained a minimal offset of  <dig>  <cit> .
significant differential expression of genes was determined by directly contrasting the control signal to the ec <dig> exposure samples and fitting a linear model through the data all using the limma package for r/bioconductor  <cit>  after which pp-values were adjusted for multiple testing according to the benjamini-hochberg method.
a gene set enrichment analysis was performed in blast2go using a fisher’s exact test again using a p-value correction with the false discovery rate  method by benjamini-hochberg  <cit> .
RESULTS
transcriptome sequencing and assembly
quality filtering of the reads was performed before assembly by applying default parameters for nontrace sequences as described by vera et al.
<cit> .
after the pyrosequencing run, 68%  of the total number of wells passed the quality control filter, representing  <dig> , <dig> bases of raw ngs data.
the reads obtained from this run had median and average lengths of  <dig> and  <dig> bases, respectively; and the longest read was  <dig>  bases long .
the obtained large number of bases into long reads indicated high quality transcriptome-wide sequence information.
the assembly was performed on all  <dig> , <dig> quality-controlled reads summing up  <dig> % of raw ngs data , after trimming of adaptor sequences .
newbler aligned  <dig> % and assembled  <dig> % of all reads, resulting in  <dig>  contigs and  <dig>  singletons .
the median coverage was  <dig>  average was  <dig> and maximum coverage was  <dig>
e. crypticus
transcriptome sequences
1after quality filtering steps and removal of outliers such as adaptor sequences and repeats.
2n <dig> is a weighted median statistic, such that 50% of all bases are contained in sequences not inferior to n <dig> length.
in the present study the largest singleton was  <dig>  bp long.
the largest contig was  <dig>  bp long and  <dig>  contigs were over  <dig>  bp long.
n <dig> is a weighted median statistic, such that 50% of all bases are contained in sequences not inferior than the n <dig> length.
the n <dig> for contig length was  <dig>  bp.
the assembly procedure using newbler considers alternative splicing, and this resulted in the integration of  <dig>  contigs into  <dig>  isotigs representing candidate transcripts.
furthermore, the assembler indicated that the isotigs are represented by  <dig>  isogroups, defined as candidate gene sequences.
the n <dig> for isotigs was  <dig>  bp and the largest isotig was  <dig>  bp long.
the largest contig count was  <dig> per single isotig and  <dig> contigs per single isogroup.
moreover 80% of isotigs and 90% of isogroups comprised one single contig.
quality control of transcriptome completeness
several quality control analyses were performed in order to assess quality and completeness of the present transcriptome.
the core eukaryotic genes mapping approach  contains a core set of  <dig> single copy genes present in all eukaryotic genomes that can be used to assess the completeness of genome and transcriptome assemblies.
the drosophila melanogaster reference set of core eukaryotic genes   <cit>  was used in blastx searches against the e. crypticus transcriptome.
at a threshold expect value of 10- <dig>   <dig> of the  <dig> cegs could be retrieved, indicating that the transcriptome is about 86% complete .
this percentage increased slightly to 92% when we used a more relaxed e-value threshold of 10- <dig>  we also asked the question whether the lack of annotation  is due to technical errors of the sequencing run or assembly errors.
therefore, we verified whether the annotation process was influenced by the quality of the sequence data.
to that end, we compared open reading frame  length, mean sequence length, percentage of start codons and mean gc content of contigs and singletons according to the quality control methods and guidelines for transcriptome sequence data described by gibson et al.
<cit> .
table  <dig> summarizes the outcome of this analysis.
firstly, annotated sequences were on average significantly longer than unannotated .
likewise, sequence length of annotated orfs  was significantly longer than unannotated orfs .
moreover, orfpredictor predicted a higher percentage of start codons in annotated sequences when compared to unannotated: as much as 84% of all contigs contained a start codon.
gc content was significantly higher in annotated sequences when compared to unannotated .
a study comparing anopheles’ and drosophila’s genomes with their corresponding proteomes showed that orphan  sequences are indeed protein-coding, and mostly code for shorter proteins when compared to currently annotated proteins  <cit> .
therefore, unannotated orfs are expected to be more abundant among shorter predicted orfs as summarized in table  <dig> and especially in relation to contigs .
functional annotation of the e. crypticus transcriptome
as this represents the first e. crypticus transcriptome sequencing study, it was crucial to include all sequence data in the analyses, to avoid losing a priori valuable biological information.
table  <dig> summarizes the functional annotation of e. crypticus’ transcriptome .
the annotation process retrieved information on gene ontology  terms, functional protein domains via interproscan  <cit> , enzyme codes and enzymatic pathway information from kyoto encyclopedia of genes and genomes .
from the total number of singletons and contigs, a blastx match was obtained for 29% and 50% respectively .
this implies that between 50%  and 70%  of the sequences do not show homology to any other sequence present in the investigated databases.
yet, over 77% of unknown contigs and 65% of unknown singletons are transcribed as determined by microarray analysis , suggesting that they may represent protein-coding genes as well as non-coding trasncripts such as micrornas and long non-coding rnas.
e. crypticus
transcriptome annotation
gene ontology terms were assigned to 25% of singletons and 44% of contigs; overall more than  <dig>  sequences were mapped.
in total  <dig>  sequences  received a full annotation .
figure  <dig> shows the blastx top-hits for both singletons and contigs, listing the top species for which highest sequence similarity was obtained against e. crypticus sequences.
orthomcl analysis  <cit>  revealed that  <dig>  contigs show orthology to one or more organisms in the orthomcl database.
within singletons we revealed orthology for  <dig>  sequences.
most singletons and contigs showed high similarity with invertebrate deuterostomes such as: the chordate branchiostoma floridae , the hemichordate acorn worm saccoglossus kowalevskii , and the purple sea urchin strongylocentrotus purpuratus  .
blastx top-hits for  <dig>  singletons and  <dig>  contigs were obtained by comparison with b. floridae, while top-hits for  <dig>  singletons and  <dig>  contigs were retrieved by comparison with s. kowalevskii.
the large number of matches between e. crypticus and chordates confirms the suggestion that the lophotrochozoan gene complement represents a rather complete bilaterian genome and that relatively few genes were lost in comparison to ecdysozoan genomes  <cit> .
this was recently confirmed by investigating the syntheny in three spiralian genomes .
they showed that the analysed worm and snail genomes are more similar to invertebrate deuterostomes when considering macrosynteny , microsynteny  of  <dig> gene clusters including the hox gene complex  <cit> , despite the fact that molluscs and annelids are evolutionary more closely related to protostomes such as arthropods, nematodes and flatworms.
additional file  <dig> illustrates the general distribution of go terms across the directed acyclic graph  structure at level  <dig>  which is based on the percentages of mapped sequences for singletons and contigs with a cutoff at 7%.
in total, thirty go terms represented the most abundant annotation of singletons and contigs.
when contigs are compared to singletons, we observe that go term distribution is balanced with regard to biological process as well as cellular compartment.
in contrast, go terms related to the molecular functions cation binding and nucleotide binding are clearly more present within contigs, while go terms related to the molecular functions protein binding, small molecule binding, hydrolase activity, ion binding, transferase activity and nucleic acid binding are more present within singletons.
this suggests that certain molecular functions may be associated with transcript length, but further research should provide more detailed information on this.
finally, the enzyme numbers, retrieved from the annotation of contigs and singletons, were submitted to ipath interactive pathway explorer to verify to what extent metabolic pathways and biosynthesis pathways are represented in the transcriptome.
as shown in additional file  <dig>  the present e. crypticus transcriptome comprises most of all essential metabolic pathways.
metabolism for terpenoids is not represented as such pathways are exclusive to plants.
microarray validation and data analysis
the agilent 4x180k microarray platform has a maximum capacity of  <dig>  spots.
after quality filtering of the contigs by earray this turned out to be sufficient to include a probe for each contig and singleton in duplicate on the array.
sequence orientation was verified by orf orientation upon translation, and showed that directional pyrosequencing was successful, because only about ~ <dig> sequences were identified to be 3&vprime; to 5&vprime; orientated.
we made sure that all annotated contigs and singletons were represented on the array.
the e. crypticus array was tested in a gene expression profiling experiment, in which we studied differential expression in between control animals and animals exposed to the zinc ec <dig> on reproduction.
this also allowed us to quantify the percentage of contigs and singletons that are actually transcribed in the control and in a stressful condition.
since the array contains  <dig> negative control spots , we quantified a mean background hybridization signal for the cy- <dig> and cy- <dig> channels separately by taking the average intensity of these spots per channel for all microarrays and calculating the overall average intensity for the red and the green channel.
subsequently, 1x the standard deviation  of this average intensity  was used as a cut-off value to verify if a sequence was actually transcribed.
this resulted in 67% of all probes showing transcriptional activity, which is comparable to previously reported percentage of detectable levels of condition-dependent transcription  <cit> .
we calculated absolute average intensity values of all contigs and compared that to average absolute intensity values of all singletons .
the average intensity values for contigs was almost  <dig> fold higher than for singletons, indicating that the majority of probes derived from singletons show lower average signal intensity.
indeed, a t-test confirmed that significantly more contigs were detected when compared to singletons .
unknown sequences may include intronic regions or untranscribed sequences that remained in the library despite dnase i treatment followed by reverse transcription and cdna synthesis.
if this is the case, probes derived from such sequences do not represent functional genes and will not show a hybridization signal.
remarkably, 77% of probes derived from unknown contigs and 66% probes derived from unknown singletons showed detectable levels of transcription.
this suggested that most unannotated sequences belong to genomic regions that are being transcribed.
additional file  <dig> shows ma plots of all probes , contigs  and singletons .
log <dig> average intensity values range between  <dig>  and  <dig>  which is very well in line with average microarray hybridization data  <cit> .
the fact that m is independent of a over a large range demonstrates sufficient linearity of the responses over a large dynamic range of microarrays signal intensity.
still, a substantial collection of probes seem to show very low average intensity values.
as mentioned above, these low intensity probes were most abundant among singletons.
exposure to zinc resulted in  <dig> significantly regulated transcripts , of which  <dig> are annotated .
figure 4b represents the volcano plot of this experiment.
among these transcripts, two were annotated as zinc efflux proteins and are both up-regulated between  <dig>  and  <dig>  fold .
additionally, a zinc transporter  <dig> was identified , which was  <dig>  fold up-regulated.
also, a zinc cadmium resistance protein  was  <dig>  fold up-regulated.
these up-regulated transcripts are associated with zinc stress, indicating that an active efflux mechanism is induced aiming to remove excess zinc ions  <cit> .
in order to retrieve biological processes affected by zinc we performed go enrichment analysis  <cit> .
table  <dig> lists the  <dig> go terms that were significantly enriched among the significantly regulated and annotated transcripts.
the relatedness of these terms is drawn in a cyclical graph in additional file  <dig>  the process of protein translation, assembly and localisation targeted to membranes is represented by  <dig> go terms.
it is tempting to speculate that zinc exposure adversely affects protein metabolism, but experimental evidence is needed to confirm this.
nevertheless, this is in agreement with the observation that zinc is an essential metal serving as a critical structural component of many zinc-finger proteins with diverse functions  <cit> .
the affected process of targeted localisation to membranes is also in line with the observation that efflux transporters are up-regulated, because these proteins have to be targeted and assembled in the cellular membrane to function properly.
1
2
3
1false discovery rate corrected significant p-value ; 2number of significantly regulated transcripts containing the particular go id; 3number of transcripts on the array containing the particular go id.
unexpectedly, five go terms representing biological processes associated with viral transcription and reproduction were significant in the enrichment analysis.
closer inspection of the transcripts from these go terms showed that most of these genes are associated with transposon activity .
jockey is a non-long terminal repeat retrotransposon , commonly found in invertebrates and amplifies via an rna intermediate  <cit> .
the autonomous retrotransposition mechanism of n-ltrs is facilitated by a single open reading frame containing polymerases with zinc-finger domains for nucleic acid binding  <cit> .
apparently, zinc exposure inactivates retrotransposition, since all these transcripts are down-regulated .
our results suggest that excess zinc ions may impair proper nucleic acid binding facilitated by the zinc-finger domains of the n-ltr open reading frame.
the final three processes among the go-enriched terms  are involved in rna catabolism.
finally, we compared the e. crypticus transcriptome to transcript sequences of the sister species e. albidus available in genbank and enchybase  <cit>  by using blastx.
a total of  <dig> e. albidus sequences clustered into  <dig> unique gene clusters, of which  <dig> showed highly significant blastx similarities to contigs in the e. crypticus transcriptome.
novais et al.
<cit>  studied the impact of cadmium and zinc on e. albidus and also included a 4 day ec <dig> zinc treatment.
although this microarray study only assessed a limited amount  of e. albidus sequences we were able to retrieve  <dig> orthologous transcripts that show highly comparable significant gene regulation upon zn exposure.
two of the ten transcripts were cysteine rich proteins, which contained metallothionein- <dig> clusters denoting putative metallothionein genes.
other genes of interest are tropomyosin- <dig> and calcium binding alpha-actinin, which are both involved in the cell-skeleton integrity.
moreover, homologs of several gene families show comparable significant regulation in e. albidus as well as e. crypticus zn exposure studies.
seven members of the ribosomal protein family are significantly up regulated in e. crypticus, while two ribosomal protein genes are up regulated in e. albidus.
also, both species show up regulation of actin, glycine n-methyltransferase and nicotinic acetylcholine receptor alpha subunit precursor homologs .
this indicates that substantial degree of conservation of transcriptional regulation in response to zn stress can be observed among related species.
CONCLUSIONS
in conclusion we showed that the current microarray platform is functional in terms of identification of transcripts associated with toxicity.
zinc-exposed animals were not included in the initial rna pool.
still, we identified differentially expressed genes responding to zinc toxicity, which were identified in previous biochemical studies to be important in zinc detoxification.
additionally, microarray data suggest that transposon activity is altered by zinc toxicity, which are highly interesting avenues for further  research.
however, this requires a more elaborate experimental design that addresses developmental programming, dose–response measurements and inter-individual variation.
also.
we deduce from the scientific literature that expression microarrays are nowadays more and more supplanted by rnaseq analysis.
according to robbens et al.
<cit>  lack of standardization and cross-platform validation prevents acceptance of genomic data integration into risk assessment.
however, more recently, the microarray quality control consortium provided such standardization and cross-validation data specifically for microarray platforms  <cit> .
this led to acceptance of this tool in clinical and pharmaceutical testing, while rnaseq technology currently lacks such high level of standardization.
it also initiated case studies that successfully integrated microarray data into environmental risk assessment  <cit> , which indicates that microarray analysis will remain an important tool, especially in 21st century risk assessment practice.
availability of supporting data
raw sequence read data were submitted to the ncbi sequence read archive under accession number .
the assembled e. crypticus contigs were submitted to the ncbi transcriptome shotgun assembly  sequence database under the tsa accession number .
raw microarray data and its miame compliant metadata were deposited at ncbi gene-expression omnibus under accession number  and platform id .
all three datasets are organized under bio-project .
competing interests
the authors declare that they have no competing interests.
authors’ contributions
mpcf: performed all exposures and rna isolations, participated in assembly and annotation, participated in microarray construction, performed data analysis, participated in design and coordination of the study, drafted the manuscript.
teb: performed microarray construction, designed microarray experiment and performed data analysis, including microarray analysis, statistical analysis and pathway analysis, participated in drafting the manuscript.
jkc: facilitated  <dig> pyrosequencing, assembly and annotation, participated in design and coordination of the study, helped drafting the manuscript.
rv: performed microarray hybridizations and scanning.
camg: participated in design of experimental treatments, helped planning the study, helped drafting the manuscript.
nms: helped planning the study, helped drafting the manuscript.
amvms: helped planning the study, helped drafting the manuscript.
mjba: participated in design and coordination of the study, helped drafting the manuscript.
dr: designed and coordinated the study, drafted the manuscript, participated in microarray construction, participated in annotation and data analysis.
all authors read and approved the final manuscript.
supplementary material
additional file 1
detailed description of all treatments, including data on rna quantity and quality.
click here for file
additional file 2
description of compounds applied in the chemical exposure with nominal sublethal effect concentrations and cas numbers.
click here for file
additional file 3
length distribution of the raw sequence reads.
click here for file
additional file 4
core eukaryotic genes mapping approach  to assess the completeness ofe.
crypticus’ transcriptome.
the  <dig> core proteins from drosophila were used for homology searches in the transcriptome.
click here for file
additional file 5
summary of annotations as performed by the transcriptome annotation pipeline at cgb, indiana university .
the file is subdivided in annotations for singletons, contigs and isotigs.
click here for file
additional file 6
distribution of most abundant go terms fore.
crypticussingletons and contigs.
level  <dig> go terms, obtained for biological processes , molecular functions  and cellular compartments , are represented by the percentage of sequences mapped  to singletons or contigs.
click here for file
additional file 7
resulting graph of ipath metabolic pathway explorer.
all enzyme codes from singleton & contig annotations were used as input.
blue-grey edges represent all enzymatic reactions in ipath explorer.
red edges represent enzymatic reactions deduced from annotation in the e. crypticus transcriptome.
blue-grey dots: chemical compound.
click here for file
additional file 8
mean average  plots for all probes , contigs  and singletons .
click here for file
additional file 9
list of significantly zinc-regulated  transcripts with annotations.
click here for file
additional file 10
go acyclic graph of the significantly zn-regulated annotated transcripts.
click here for file
acknowledgements
thanks to dr. k. mokaitis, dr j.-h. choi, z. smith, j. ford, a. buechlein at the centre for genomics and bioinformatics, indiana university, for sequencing and bioinformatics; janine mariën and rudo a. verweij at vu university amsterdam, for technical support in the lab.
this study was financially supported by feder funding through compete programa operacional factores de competitividade, and by portuguese funding through fundação para a ciência e tecnologia , within the research project nanoka fcomp-01-0124-feder- <dig> , and through the fct phd grant to m.p.
castro-ferreira .
the ecogenomics consortium “ecolinc” research project was funded by netherlands organisation for scientific research .
additional support was provided by the metacyt initiative of indiana university, funded in part through a major grant from the lilly endowment, inc.

 the ecologically relevant soil-living worm enchytraeus crypticus has been widely used as model species in soil ecotoxicology for over two decades, both under international organization for standardization and organization for economic co-operation and development guidelines. ecotoxicological tests generally focus on survival and reproduction endpoints, i.e. organism effects, whereas in ecotoxicogenomics, responses at the transcriptome level are studied aiming to reconstitute the mechanisms preceding effects at organismal and/or population level. in-depth information about the molecular mechanisms is provided allowing a better understanding of molecular pathways activated by toxicity, which will not be possible when using the classical tests alone. previously, transcript information with classical sanger sequencing was generated for the test species enchytraeus albidus and has been very useful to further explore the response mechanisms of enchytraeids to various stressors, such as pesticides, heavy metals, nanoparticles, and soil properties. however, the drawback of the current e. albidus microarray is that it contains a low number of transcripts, so a full transcriptomic assessment of stress response cannot be executed for this model. thanks to the progress in sequencing technology, transcriptome and genome information can now be obtained for any non-genomic model species. transcriptome sequencing is of special interest, since it is directly targeted to obtain sequence information from transcribed regions in the genome. to pursue a comprehensive transcriptome of an enchytraeid test species we selected e. crypticus instead of the more widely used e. albidus, despite the seemingly lower sensitivity to toxicants of this species. several advantages over the use of e. albidus have been recently identified, such as easy culturing conditions in agar media, higher reproductive rate, lack of cryptic speciation, and a shorter reproduction cycle. given the clear advantages of using e. crypticus for soil quality assessment, the study presented here aims at studying its transcriptome by de novo transcriptome pyrosequencing and annotation. a previous transcriptome study showed that it is beneficial to pool rna from multiple treatments in order to obtain a wider range of expressed transcript sequences. thus, we decided to set up a similar approach to get transcriptome-wide information from e. crypticus and sampled mrna from multiple exposure conditions, including chemical, physicochemical and developmental stages. the quality of the assembled transcriptome was critically assessed according to the guideline provided by gibson et al. by evaluating sequence quality, assembly quality, overall transcriptome completeness, and transcriptional activity. subsequently, a high-density microarray platform was constructed and validated using zinc-exposed organisms. this platform will be a valuable resource to conduct gene expression studies in future ecotoxicogenomic research. methods test organism culturing conditions and experiment overview the annelid enchytraeus crypticus, commonly known as potworm, has been cultured for several years in laboratory at the vu university amsterdam. cultures are kept in agar substrate prepared with aqueous soil extract, fed ad libitum with ground oatmeal, at 16°c with 75% relative humidity and 16/8 h light/dark photoperiod. for this study we selected, for each biological replicate per exposure, 5– enchytraeus crypticus adults with developed clitella, similar size, and an average fresh weight of mg per adult. in total, different experimental conditions were included: chemical exposures temperature treatments, soil moisture treatments, ph treatments, different developmental stages, plus associated control conditions. test soils the reference natural lufa soil was used in all tests, except for the ph treatments, which were performed using oecd artificial soil. lufa soil is a loamy sand soil sieved at 2 mm with ph, cation exchange capacity 10 meq/100 g, maximum water-holding capacity of 46%, organic matter content of 4%, grain size distribution of 6% clay, 14% silt and 80% sand. for testing, control moisture was adjusted to 50% of the whcmax. the oecd artificial soil was prepared according to the standard procedure and further used to test various phs adjusting it by adding predetermined amounts of calcium carbonate. treatments toxicity tests were performed following the standard guidelines. chemical test concentrations were based on previous results from which the 3-week reproduction ec was estimated. this ec was used as a reference exposure concentration in the short-term experiments from which the mrna pool was derived. test chemicals included cadmium chloride hemipentahydrate, carbendazim, phenanthrene, pentachloroaniline and -dichloroaniline. chemical formulas, chemical abstract service numbers and concentrations are described in additional file non-water soluble chemicals were solubilized in acetone. for the acetone control and acetone chemical spiking the following was performed: gr of dry lufa soil was spiked with acetone or the appropriate organic-acetone solution for phe, car, pca and dca, equilibrated for 24 h, left to evaporate overnight, mixed with dry soil and homogeneously moistened to 50% of the whcmax. the soil batches were used immediately after being prepared. for the metal exposures, dry lufa soil was spiked with the specific metal-aqueous solution to reach the nominal concentration and 50% of the whcmax in the soil batch. the soil was kept under standard conditions while equilibrating for 21 days prior to exposure. two replicates were prepared per test condition. per replicate, twenty adult e. crypticus were introduced into each glass vial containing 30 g of soil prepared as described above. food was supplied and vessels were closed with perforated aluminum foil. the exposures lasted 48 h at 20°c, with 75% relative humidity and 16/8 h light/dark photoperiod. after 48 h of exposure, the individuals in each replicate were sampled, carefully rinsed in deionised water, transferred to a cryotube, snap-frozen with liquid nitrogen and stored at -80°c. zinc exposures were performed in order to validate the microarray platform. the ec for inhibition of reproduction was used. for the microarray analysis, e. crypticus adults were exposed to 145 mg zn/kg soil and control soil for 4 days using four replicates, containing 5– individuals per replicate. temperature, moisture and ph treatments were performed according to the guidelines, with standard test conditions as described in the previous paragraph, except for the specific factor being studied. the temperature exposures lasted, or 96 h, at the temperatures - and 40°c. the moisture exposures lasted or 96 h, using soils moistened to and 90% of whcmax. the ph exposures lasted 48 h, using soils with ph  and enchytraeus crypticus samples in distinct developmental stages were taken from synchronization of cultures, which were obtained by transferring thirty adults with developed clitella to fresh culture boxes that were kept at identical culturing conditions as the source culture. every 6 hours, cocoons were collected, immediately snap-frozen with liquid nitrogen and stored at -80°c. twenty synchronized cocoons were carefully transferred into fresh culture media and kept at 20°c. after twenty days twenty juveniles were collected into a cryotube, immediately snap-frozen with liquid nitrogen and stored at -80°c. rna isolation, rna pool and cdna library construction total rna was collected from the individual samples of enchytraeus crypticus as listed in additional file all rna extractions were performed using the sv total rna isolation system following manufacturer’s protocol. all samples were homogenized with liquid nitrogen in the cryotube before rna extraction. addition of sterilized quartz sand to the liquid nitrogen was needed to disrupt and homogenize cocoons sufficiently. the overall quality of rna samples was evaluated in terms of purity and integrity of rna by means of a nanodrop nd- uv–vis spectrometer and agarose gel electrophoresis. rna sample quality was verified regarding high rna concentration, absorbance ratios a260/a in the range -, and a260/a above. samples with lower absorbance ratio were ethanol-precipitated in order to improve the quality. equivalent amounts of rna mass per test condition were pooled together, with a total of 10 μg rna from all samples of e. crypticus. figure  illustrates the proportion of rna input per test condition represented in the pool that was used for sequencing. normalized 454-sequencing libraries were constructed from an equal-molar pool of rna obtained from the unique exposure samples described above using the procedures described in meyer et al.. after the final purification step the library was stored at -20°c until sequencing. transcriptome sequencing and assembly cytochrome oxidase i and 18s rrna gene sequencing of individuals yielded only one polymorphic site, providing strong evidence that this strain is highly inbred. thus, the assembly procedure is not hampered by the presence of genetic variation within the cdna library. this library was submitted to one full-plate sequencing run using a roche gs flx pyrosequencing instrument using titanium chemistry, following manufacturer’s protocol and methods previously described. after sequencing, the generated sequence reads were assembled using newbler in de novo mode and default parameters. transcriptome annotation transcriptome annotation was performed through the isga transcriptome analysis pipeline, available from the center for genomics and bioinformatics at indiana university. first, sequence homology to known metazoan proteins was obtained by submitting contigs and singletons to blastx searches against ncbi’s non-redundant database and dbest. moreover, protein domains were identified among the frame translations of the assembled sequences using pfam, tigrfam and hmmer searches. open reading frames were determined with the orfpredictor software on the proteomics server of youngstown state university. orthology and paralogy was assigned by blastx against orthomcl databases. assembly quality was tested by retrieving blastx hits against the drosophila orthologs in the cegma core eukaryotic genes dataset. transcriptome quality assessment was performed according to the guidelines described by gibson et al.. in short, open reading frame lengths, est nucleotide lengths and gc contents of annotated and unannoted ests were compared using a welch’s t-test in r version. furthermore, transcript abundance was compared between contigs and singletons using welch’s t-test. finally, the interactive pathway explorer ipath was used to identify metabolic, regulatory and biosynthesis pathways based on enzyme code numbers. probe design for microarray construction oligonucleotide probes were designed according to the thermodynamics/technical features of gene expression microarrays within the earray webtool of agilent. the workflow is summarized into these consecutive steps: the target sequences were selected by length threshold ≥ bases, converted into sense orientation, and then submitted to earray for probe design. in order to avoid cross-hybridization, the probes that could cross-hybridize with any sequence other than the specific target, were discarded. a final microarray design was created from the probe group that were printed in duplicate onto a x k platform. the array design was supplemented with additional agilent control probes, so that in total probes were present. microarray validation and data analysis rna samples were extracted from control and ec zn-exposed e. crypticus adults. the rna was then precipitated with ethanol to remove any guanidine thiocyanate traces from the sv rna extraction procedure and to increase its concentration. rna concentration and quality was determined on the nanodrop nd spectrophotometer and bioanalyzer. rna samples were labeled using the low input quick amp kit according to the manufacturer’s instructions and as described. four replicates were labeled using cy-3/cy- dye swaps. labeled samples were hybridized to the microarray as described above. microarrays were scanned using an agilent microarray scanner and raw data was extracted with feature extraction software version.. microarray raw intensities were loess normalized with a background subtraction according to the edwards method which contained a minimal offset of. significant differential expression of genes was determined by directly contrasting the control signal to the ec exposure samples and fitting a linear model through the data all using the limma package for r/bioconductor after which pp-values were adjusted for multiple testing according to the benjamini-hochberg method. a gene set enrichment analysis was performed in blast2go using a fisher’s exact test again using a p-value correction with the false discovery rate method by benjamini-hochberg.  transcriptome sequencing and assembly quality filtering of the reads was performed before assembly by applying default parameters for nontrace sequences as described by vera et al.. after the pyrosequencing run, 68% of the total number of wells passed the quality control filter, representing, bases of raw ngs data. the reads obtained from this run had median and average lengths of and bases, respectively; and the longest read was bases long. the obtained large number of bases into long reads indicated high quality transcriptome-wide sequence information. the assembly was performed on all, quality-controlled reads summing up % of raw ngs data, after trimming of adaptor sequences. newbler aligned % and assembled % of all reads, resulting in contigs and singletons. the median coverage was average was and maximum coverage was e. crypticus transcriptome sequences 1after quality filtering steps and removal of outliers such as adaptor sequences and repeats. 2n is a weighted median statistic, such that 50% of all bases are contained in sequences not inferior to n length. in the present study the largest singleton was  bp long. the largest contig was  bp long and contigs were over  bp long. n is a weighted median statistic, such that 50% of all bases are contained in sequences not inferior than the n length. the n for contig length was  bp. the assembly procedure using newbler considers alternative splicing, and this resulted in the integration of contigs into isotigs representing candidate transcripts. furthermore, the assembler indicated that the isotigs are represented by isogroups, defined as candidate gene sequences. the n for isotigs was  bp and the largest isotig was  bp long. the largest contig count was per single isotig and contigs per single isogroup. moreover 80% of isotigs and 90% of isogroups comprised one single contig. quality control of transcriptome completeness several quality control analyses were performed in order to assess quality and completeness of the present transcriptome. the core eukaryotic genes mapping approach contains a core set of single copy genes present in all eukaryotic genomes that can be used to assess the completeness of genome and transcriptome assemblies. the drosophila melanogaster reference set of core eukaryotic genes was used in blastx searches against the e. crypticus transcriptome. at a threshold expect value of 10- of the cegs could be retrieved, indicating that the transcriptome is about 86% complete. this percentage increased slightly to 92% when we used a more relaxed e-value threshold of 10- we also asked the question whether the lack of annotation is due to technical errors of the sequencing run or assembly errors. therefore, we verified whether the annotation process was influenced by the quality of the sequence data. to that end, we compared open reading frame length, mean sequence length, percentage of start codons and mean gc content of contigs and singletons according to the quality control methods and guidelines for transcriptome sequence data described by gibson et al.. table  summarizes the outcome of this analysis. firstly, annotated sequences were on average significantly longer than unannotated. likewise, sequence length of annotated orfs was significantly longer than unannotated orfs. moreover, orfpredictor predicted a higher percentage of start codons in annotated sequences when compared to unannotated: as much as 84% of all contigs contained a start codon. gc content was significantly higher in annotated sequences when compared to unannotated. a study comparing anopheles’ and drosophila’s genomes with their corresponding proteomes showed that orphan sequences are indeed protein-coding, and mostly code for shorter proteins when compared to currently annotated proteins. therefore, unannotated orfs are expected to be more abundant among shorter predicted orfs as summarized in table  and especially in relation to contigs. functional annotation of the e. crypticus transcriptome as this represents the first e. crypticus transcriptome sequencing study, it was crucial to include all sequence data in the analyses, to avoid losing a priori valuable biological information. table  summarizes the functional annotation of e. crypticus’ transcriptome. the annotation process retrieved information on gene ontology terms, functional protein domains via interproscan, enzyme codes and enzymatic pathway information from kyoto encyclopedia of genes and genomes. from the total number of singletons and contigs, a blastx match was obtained for 29% and 50% respectively. this implies that between 50% and 70% of the sequences do not show homology to any other sequence present in the investigated databases. yet, over 77% of unknown contigs and 65% of unknown singletons are transcribed as determined by microarray analysis, suggesting that they may represent protein-coding genes as well as non-coding trasncripts such as micrornas and long non-coding rnas. e. crypticus transcriptome annotation gene ontology terms were assigned to 25% of singletons and 44% of contigs; overall more than sequences were mapped. in total sequences received a full annotation. figure  shows the blastx top-hits for both singletons and contigs, listing the top species for which highest sequence similarity was obtained against e. crypticus sequences. orthomcl analysis revealed that contigs show orthology to one or more organisms in the orthomcl database. within singletons we revealed orthology for sequences. most singletons and contigs showed high similarity with invertebrate deuterostomes such as: the chordate branchiostoma floridae, the hemichordate acorn worm saccoglossus kowalevskii, and the purple sea urchin strongylocentrotus purpuratus. blastx top-hits for singletons and contigs were obtained by comparison with b. floridae, while top-hits for singletons and contigs were retrieved by comparison with s. kowalevskii. the large number of matches between e. crypticus and chordates confirms the suggestion that the lophotrochozoan gene complement represents a rather complete bilaterian genome and that relatively few genes were lost in comparison to ecdysozoan genomes. this was recently confirmed by investigating the syntheny in three spiralian genomes. they showed that the analysed worm and snail genomes are more similar to invertebrate deuterostomes when considering macrosynteny, microsynteny of gene clusters including the hox gene complex, despite the fact that molluscs and annelids are evolutionary more closely related to protostomes such as arthropods, nematodes and flatworms. additional file illustrates the general distribution of go terms across the directed acyclic graph structure at level which is based on the percentages of mapped sequences for singletons and contigs with a cutoff at 7%. in total, thirty go terms represented the most abundant annotation of singletons and contigs. when contigs are compared to singletons, we observe that go term distribution is balanced with regard to biological process as well as cellular compartment. in contrast, go terms related to the molecular functions cation binding and nucleotide binding are clearly more present within contigs, while go terms related to the molecular functions protein binding, small molecule binding, hydrolase activity, ion binding, transferase activity and nucleic acid binding are more present within singletons. this suggests that certain molecular functions may be associated with transcript length, but further research should provide more detailed information on this. finally, the enzyme numbers, retrieved from the annotation of contigs and singletons, were submitted to ipath interactive pathway explorer to verify to what extent metabolic pathways and biosynthesis pathways are represented in the transcriptome. as shown in additional file the present e. crypticus transcriptome comprises most of all essential metabolic pathways. metabolism for terpenoids is not represented as such pathways are exclusive to plants. microarray validation and data analysis the agilent 4x180k microarray platform has a maximum capacity of spots. after quality filtering of the contigs by earray this turned out to be sufficient to include a probe for each contig and singleton in duplicate on the array. sequence orientation was verified by orf orientation upon translation, and showed that directional pyrosequencing was successful, because only about ~ sequences were identified to be 3&vprime; to 5&vprime; orientated. we made sure that all annotated contigs and singletons were represented on the array. the e. crypticus array was tested in a gene expression profiling experiment, in which we studied differential expression in between control animals and animals exposed to the zinc ec on reproduction. this also allowed us to quantify the percentage of contigs and singletons that are actually transcribed in the control and in a stressful condition. since the array contains negative control spots, we quantified a mean background hybridization signal for the cy- and cy- channels separately by taking the average intensity of these spots per channel for all microarrays and calculating the overall average intensity for the red and the green channel. subsequently, 1x the standard deviation of this average intensity was used as a cut-off value to verify if a sequence was actually transcribed. this resulted in 67% of all probes showing transcriptional activity, which is comparable to previously reported percentage of detectable levels of condition-dependent transcription. we calculated absolute average intensity values of all contigs and compared that to average absolute intensity values of all singletons. the average intensity values for contigs was almost fold higher than for singletons, indicating that the majority of probes derived from singletons show lower average signal intensity. indeed, a t-test confirmed that significantly more contigs were detected when compared to singletons. unknown sequences may include intronic regions or untranscribed sequences that remained in the library despite dnase i treatment followed by reverse transcription and cdna synthesis. if this is the case, probes derived from such sequences do not represent functional genes and will not show a hybridization signal. remarkably, 77% of probes derived from unknown contigs and 66% probes derived from unknown singletons showed detectable levels of transcription. this suggested that most unannotated sequences belong to genomic regions that are being transcribed. additional file shows ma plots of all probes, contigs and singletons. log average intensity values range between and which is very well in line with average microarray hybridization data. the fact that m is independent of a over a large range demonstrates sufficient linearity of the responses over a large dynamic range of microarrays signal intensity. still, a substantial collection of probes seem to show very low average intensity values. as mentioned above, these low intensity probes were most abundant among singletons. exposure to zinc resulted in significantly regulated transcripts, of which are annotated. figure 4b represents the volcano plot of this experiment. among these transcripts, two were annotated as zinc efflux proteins and are both up-regulated between and fold. additionally, a zinc transporter was identified, which was fold up-regulated. also, a zinc cadmium resistance protein was fold up-regulated. these up-regulated transcripts are associated with zinc stress, indicating that an active efflux mechanism is induced aiming to remove excess zinc ions. in order to retrieve biological processes affected by zinc we performed go enrichment analysis. table  lists the go terms that were significantly enriched among the significantly regulated and annotated transcripts. the relatedness of these terms is drawn in a cyclical graph in additional file the process of protein translation, assembly and localisation targeted to membranes is represented by go terms. it is tempting to speculate that zinc exposure adversely affects protein metabolism, but experimental evidence is needed to confirm this. nevertheless, this is in agreement with the observation that zinc is an essential metal serving as a critical structural component of many zinc-finger proteins with diverse functions. the affected process of targeted localisation to membranes is also in line with the observation that efflux transporters are up-regulated, because these proteins have to be targeted and assembled in the cellular membrane to function properly. 1 2 3 1false discovery rate corrected significant p-value; 2number of significantly regulated transcripts containing the particular go id; 3number of transcripts on the array containing the particular go id. unexpectedly, five go terms representing biological processes associated with viral transcription and reproduction were significant in the enrichment analysis. closer inspection of the transcripts from these go terms showed that most of these genes are associated with transposon activity. jockey is a non-long terminal repeat retrotransposon, commonly found in invertebrates and amplifies via an rna intermediate. the autonomous retrotransposition mechanism of n-ltrs is facilitated by a single open reading frame containing polymerases with zinc-finger domains for nucleic acid binding. apparently, zinc exposure inactivates retrotransposition, since all these transcripts are down-regulated. our results suggest that excess zinc ions may impair proper nucleic acid binding facilitated by the zinc-finger domains of the n-ltr open reading frame. the final three processes among the go-enriched terms are involved in rna catabolism. finally, we compared the e. crypticus transcriptome to transcript sequences of the sister species e. albidus available in genbank and enchybase by using blastx. a total of e. albidus sequences clustered into unique gene clusters, of which showed highly significant blastx similarities to contigs in the e. crypticus transcriptome. novais et al. studied the impact of cadmium and zinc on e. albidus and also included a 4 day ec zinc treatment. although this microarray study only assessed a limited amount of e. albidus sequences we were able to retrieve orthologous transcripts that show highly comparable significant gene regulation upon zn exposure. two of the ten transcripts were cysteine rich proteins, which contained metallothionein- clusters denoting putative metallothionein genes. other genes of interest are tropomyosin- and calcium binding alpha-actinin, which are both involved in the cell-skeleton integrity. moreover, homologs of several gene families show comparable significant regulation in e. albidus as well as e. crypticus zn exposure studies. seven members of the ribosomal protein family are significantly up regulated in e. crypticus, while two ribosomal protein genes are up regulated in e. albidus. also, both species show up regulation of actin, glycine n-methyltransferase and nicotinic acetylcholine receptor alpha subunit precursor homologs. this indicates that substantial degree of conservation of transcriptional regulation in response to zn stress can be observed among related species.  in conclusion we showed that the current microarray platform is functional in terms of identification of transcripts associated with toxicity. zinc-exposed animals were not included in the initial rna pool. still, we identified differentially expressed genes responding to zinc toxicity, which were identified in previous biochemical studies to be important in zinc detoxification. additionally, microarray data suggest that transposon activity is altered by zinc toxicity, which are highly interesting avenues for further research. however, this requires a more elaborate experimental design that addresses developmental programming, dose–response measurements and inter-individual variation. also. we deduce from the scientific literature that expression microarrays are nowadays more and more supplanted by rnaseq analysis. according to robbens et al. lack of standardization and cross-platform validation prevents acceptance of genomic data integration into risk assessment. however, more recently, the microarray quality control consortium provided such standardization and cross-validation data specifically for microarray platforms. this led to acceptance of this tool in clinical and pharmaceutical testing, while rnaseq technology currently lacks such high level of standardization. it also initiated case studies that successfully integrated microarray data into environmental risk assessment, which indicates that microarray analysis will remain an important tool, especially in 21st century risk assessment practice. availability of supporting data raw sequence read data were submitted to the ncbi sequence read archive under accession number. the assembled e. crypticus contigs were submitted to the ncbi transcriptome shotgun assembly sequence database under the tsa accession number. raw microarray data and its miame compliant metadata were deposited at ncbi gene-expression omnibus under accession number and platform id. all three datasets are organized under bio-project. competing interests the authors declare that they have no competing interests. authors’ contributions mpcf: performed all exposures and rna isolations, participated in assembly and annotation, participated in microarray construction, performed data analysis, participated in design and coordination of the study, drafted the manuscript. teb: performed microarray construction, designed microarray experiment and performed data analysis, including microarray analysis, statistical analysis and pathway analysis, participated in drafting the manuscript. jkc: facilitated pyrosequencing, assembly and annotation, participated in design and coordination of the study, helped drafting the manuscript. rv: performed microarray hybridizations and scanning. camg: participated in design of experimental treatments, helped planning the study, helped drafting the manuscript. nms: helped planning the study, helped drafting the manuscript. amvms: helped planning the study, helped drafting the manuscript. mjba: participated in design and coordination of the study, helped drafting the manuscript. dr: designed and coordinated the study, drafted the manuscript, participated in microarray construction, participated in annotation and data analysis. all authors read and approved the final manuscript. supplementary material additional file 1 detailed description of all treatments, including data on rna quantity and quality. click here for file additional file 2 description of compounds applied in the chemical exposure with nominal sublethal effect concentrations and cas numbers. click here for file additional file 3 length distribution of the raw sequence reads. click here for file additional file 4 core eukaryotic genes mapping approach to assess the completeness ofe. crypticus’ transcriptome. the core proteins from drosophila were used for homology searches in the transcriptome. click here for file additional file 5 summary of annotations as performed by the transcriptome annotation pipeline at cgb, indiana university. the file is subdivided in annotations for singletons, contigs and isotigs. click here for file additional file 6 distribution of most abundant go terms fore. crypticussingletons and contigs. level go terms, obtained for biological processes, molecular functions and cellular compartments, are represented by the percentage of sequences mapped to singletons or contigs. click here for file additional file 7 resulting graph of ipath metabolic pathway explorer. all enzyme codes from singleton & contig annotations were used as input. blue-grey edges represent all enzymatic reactions in ipath explorer. red edges represent enzymatic reactions deduced from annotation in the e. crypticus transcriptome. blue-grey dots: chemical compound. click here for file additional file 8 mean average plots for all probes, contigs and singletons. click here for file additional file 9 list of significantly zinc-regulated transcripts with annotations. click here for file additional file 10 go acyclic graph of the significantly zn-regulated annotated transcripts. click here for file acknowledgements thanks to dr. k. mokaitis, dr j.-h. choi, z. smith, j. ford, a. buechlein at the centre for genomics and bioinformatics, indiana university, for sequencing and bioinformatics; janine mariën and rudo a. verweij at vu university amsterdam, for technical support in the lab. this study was financially supported by feder funding through compete programa operacional factores de competitividade, and by portuguese funding through fundação para a ciência e tecnologia, within the research project nanoka fcomp-01-0124-feder-, and through the fct phd grant to m.p. castro-ferreira. the ecogenomics consortium “ecolinc” research project was funded by netherlands organisation for scientific research. additional support was provided by the metacyt initiative of indiana university, funded in part through a major grant from the lilly endowment, inc.
an initial investigation of the e. crypticus transcriptome including an associated microarray platform for future studies proves to be a valuable resource to investigate functional genomics mechanisms of toxicity in soil environments and to annotate a potentially large number of lineage specific genes that are responsive to environmental stress conditions. concerning assembly quality, 84% of the contig sequences contain an open reading frame with a start codon while e. crypticus homologs were identified for 92% of the core eukaryotic genes. unexpectedly, the microarray data suggest that zinc exposure alters retro transposon activity in the e. crypticus genome. in total, significantly differentially regulated transcripts were identified upon zinc exposure. several significantly regulated genes exerted predicted functions associated with zinc stress. overall, 70% of all probes signaled expression above background levels. an agilent 180 k microarray platform was designed and validated by hybridizing cdna from 4 day zinc- exposed e. crypticus to the concentration corresponding to 50% reduction in reproduction after three weeks. the soil worm enchytraeus crypticus is an ecotoxicology model species that, until now, was without genome or transcriptome sequence information. the present research aims at studying the transcriptome of enchytraeus crypticus, sampled from multiple test conditions, and the construction of a high-density microarray for functional genomic studies.
 the soil worm enchytraeus crypticus is an ecotoxicology model species that, until now, was without genome or transcriptome sequence information. the present research aims at studying the transcriptome of enchytraeus crypticus, sampled from multiple test conditions, and the construction of a high-density microarray for functional genomic studies.  over million cdna sequence reads were obtained representing million nucleotides. after assembly, contigs and singletons were obtained, from which 44% and 25% are annotated as protein-coding genes, respectively, sharing homology with other animal proteomes. concerning assembly quality, 84% of the contig sequences contain an open reading frame with a start codon while e. crypticus homologs were identified for 92% of the core eukaryotic genes. moreover, 65% and 77% of the singletons and contigs without known homologs, respectively, were shown to be transcribed in an independent microarray experiment. an agilent 180 k microarray platform was designed and validated by hybridizing cdna from 4 day zinc- exposed e. crypticus to the concentration corresponding to 50% reduction in reproduction after three weeks. overall, 70% of all probes signaled expression above background levels. more specifically, the probes derived from contigs showed a wider range of average intensities when compared to probes derived from singletons. in total, significantly differentially regulated transcripts were identified upon zinc exposure. several significantly regulated genes exerted predicted functions associated with zinc stress. unexpectedly, the microarray data suggest that zinc exposure alters retro transposon activity in the e. crypticus genome.  an initial investigation of the e. crypticus transcriptome including an associated microarray platform for future studies proves to be a valuable resource to investigate functional genomics mechanisms of toxicity in soil environments and to annotate a potentially large number of lineage specific genes that are responsive to environmental stress conditions. ecotoxicogenomicsnext-generation pyrosequencinginvertebratezincannelid sequencing 