MARGEN-00327; No of Pages 5 Marine Genomics xxx (2015) xxx–xxx
Contents lists available at ScienceDirect
Marine Genomics journal homepage: www.elsevier.com/locate/margen
Genomics/technical resources
Spatio-temporal transcriptome analysis in the marine snail Tegula atra along central-northern Chile (28–31°S) Fabian J. Tapia a,c, Cristian Gallardo-Escárate b,c,⁎ a b c
Department of Oceanography & COPAS South-Austral, University of Concepción, Concepción, Chile Department of Oceanography, Laboratory of Biotechnology and Aquatic Genomics, University of Concepción, Concepción, Chile Interdisciplinary Center for Aquaculture Research (INCAR), University of Concepción, Concepción, Chile
a r t i c l e
i n f o
Article history: Received 15 March 2015 Received in revised form 7 May 2015 Accepted 7 May 2015 Available online xxxx Keywords: RNA-Seq Tegula atra Snail Transcriptome Sequencing by synthesis
a b s t r a c t This study describes the results from transcriptomes sequenced by illumina technology from four populations in the marine snail Tegula atra along central-northern Chile (28–31°S) during summer and winter 2014. In silico differential expression of transcripts annotated to known proteins revealed several local patterns associated to the environmental thermal variability. Herein, northern populations evidenced lower number of genes highly regulated, while southern populations displayed opposite patterns. This relationship could suggest that northern snail populations are more adapted to high temperature variations, enabling specific genes (e.g. HSPs and antioxidant system) to maintain high transcriptional activity under controlled physiological conditions. This transcriptome response or “frontloading” strategy can significantly increase the speed of response to thermal stress, and also be a relevant molecular underpinning to explain the genomic diversity along the Chilean coast. © 2015 Elsevier B.V. All rights reserved.
1. Introduction Understanding how environmental interactions affect the patterns observed in the intertidal community is key factor to comprehend its dynamics and modulation along the Chilean coast. For instance, in many rocky intertidal habitats, the high diversity of sessile organisms is allowed by the habitat availability and biomass for foraging activities. However, invertebrate populations are strongly shaped by environmental stress (e.g. temperature, oxygen, and nutrients), that in conjunction with the larval dispersal potential creates biogeographic breaks (Haye et al., 2014; Nunez-Acuna et al., 2012). The genus Tegula is well known as back turban snail and is of the most abundant mollusks of the Californian and Chilean coasts (Valdovinos, 1999; Sagarin and Gaines, 2002). Tegula atra is distributed along the Southeast Pacific from 9°S to 44°S and it has been considered as having a low potential larval dispersal between 4 and 7 days (Kulikova and Omel'yanenko, 2000). This study describes a transcriptomic-wide expression analysis in four population of the marine snail T. atra along central-northern Chile (28–31°S) during the summer and winter 2014. This transcriptome
⁎ Corresponding author at: Department of Oceanography, Laboratory of Biotechnology and Aquatic Genomics, University of Concepción, Concepción, Chile. E-mail address:
[email protected] (C. Gallardo-Escárate).
information will provide the basis to evaluate genome local adaptation in a marine species with low dispersal potential, integrating the physical variability coupled to spatio-temporal patterns. 2. Data description 2.1. RNA isolation and sequencing To overview major transcriptomic differences and thoroughly characterize the transcriptome of T. atra, sixteen cDNA libraries were constructed from RNA pools extracted from individuals sampled from four snail populations. For each population, two biological replicates (n = 10) were evaluated. The localities sampled were Huasco (28°24S– 71°11 W), Temblador (29°28S–71°20 W), Guanaqueros (31°11S– 71°28 W) and Pta. Talca (30°54S–71°51 W) in the central-northern Chile. Total RNA was extracted from gills using the RiboPure™ Kit (Ambion®, Life Technologies™, USA) following the manufacturer's instructions. Quantity, purity, and quality of isolated RNA were measured in the TapeStation 2200 (Agilent Technologies Inc., Santa Clara, CA, USA) using the R6K Reagent Kit according to the manufacturer's instructions; samples with RIN over 8.0 were used for library preparation. Subsequently, double-stranded cDNA libraries were constructed using the TruSeq RNA Sample Preparation Kit v2 (Illumina®, San Diego, CA, USA). The biological replicates were sequenced by the MiSeq (Illumina®) platform using sequenced runs of 2 x 251 paired-end
http://dx.doi.org/10.1016/j.margen.2015.05.005 1874-7787/© 2015 Elsevier B.V. All rights reserved.
Please cite this article as: Tapia, F.J., Gallardo-Escárate, C., Spatio-temporal transcriptome analysis in the marine snail Tegula atra along centralnorthern Chile (28–31°S), Mar. Genomics (2015), http://dx.doi.org/10.1016/j.margen.2015.05.005
2
F.J. Tapia, C. Gallardo-Escárate / Marine Genomics xxx (2015) xxx–xxx
reads at Interdisciplinary Center for Aquaculture Research (INCAR), University of Concepción, Chile. The cleaned short read sequences were deposited in the Sequence Read Archive (SRA) (http://www. ncbi.nlm.nih.gov/sra/SRX958768[accn]) under the accession number SRX958768. 2.2. De novo assembly and gene annotation Nine different de novo assemblies and subsequent data analyses were performed, one for each population and summer/winter, and a general assembly considering all datasets. The overlap settings for this assembly were a mismatch cost = 2, an insert cost = 3, a minimum contig length = 200 bp, a similarity = 0.8. Finally, 208,425 contigs were obtained with an average length of 501 bp. Among these, 46,639 putative transcripts were annotated to known proteins. All contigs were annotated using Blastp against the NCBI non-redundant (nr) protein database, considering a cutoff E-value of 1E-5 (Table 1S). In addition, consensus sequences from the T. atra transcriptome were annotated according to Gene Ontology (GO) terms with the Blast2GO program (Conesa et al., 2005) (Fig. 2S). For more details of the bioinformatics analyses see Supplementary Electronic Material (S1). 2.3. Transcriptome profiles from T. atra populations RNA-Seq analysis was conducted by individually mapping the Huasco, Temblador, Guanaqueros and Pta. Talca datasets against the 208,425 generated contigs (reference transcriptome), and the expression level of each transcript was quantified in reads per kilobase of the
Table 1 Differential transcription expression (winter/summer) from T. atra populations. Feature ID Huasco Contig0082710 Contig0145892 Contig0129496 Contig0017534 Contig0011841 Contig0050105 Contig0129497 Contig0013151 Contig0050143 Contig0055098 Contig0131507 Contig0085929 Contig0009090 Contig0074059 Contig0041752 Contig0040126 Contig0005967 Contig0012610 Contig0002001 Contig0017640 Contig0056160 Contig0001606 Contig0003685 Contig0006086 Contig0006069 Contig0008403 Temblador Contig0000560 Contig0200536 Contig0027852 Contig0124160 Contig0050143 Contig0027269 Contig0048336 Contig0003101 Contig0005975 Contig0037842 Contig0065331 Contig0006086 Contig0055019 Contig0098136 Contig0046761 Contig0001089 Contig0003685 Contig0028838 Contig0037487 Contig0003405 Contig0017640
Fold–change
87,09 81,92 74,17 64,95 52,03 50,93 49,03 39,85 35,61 21,77 21,40 17,57 –21,68 –26,76 –27,10 –27,32 –32,52 –38,62 –49,05 –65,71 –82,65 –87,39 –88,07 –94,84 –140,91 –387,51 182,6 63,2 38,7 38,7 32,9 20,6 20,1 19,6 –18,4 –20,8 –20,8 –26,0 –27,1 –29,1 –38,8 –47,5 –66,0 –87,2 –182,5 –189,0 –353,7
Description
Perlucin [Crassostrea gigas] Golgi–associated plant pathogenesis–related protein 1 [Crassostrea gigas] X–box binding protein [Haliotis discus discus] RecName: Full=Glycine, alanine and asparagine–rich protein; Flags: Precursor Putative ferric–chelate reductase 1–like protein [Crassostrea gigas] Invertebrate–type lysozyme [Haliotis discus discus] X–box binding protein [Haliotis discus discus] Mitochondrial malate dehydrogenase precursor, partial [Osilinus lineatus] Invertebrate–type lysozyme [Haliotis discus discus] High mobility group protein B3 [Crassostrea gigas] L–ascorbate oxidase [Crassostrea gigas] Perlucin [Crassostrea gigas] Inhibitor of apoptosis protein [Crassostrea gigas] Polyubiquitin, partial [Solen grandis] Ankyrin–2 [Crassostrea gigas] IkB–like protein [Haliotis discus discus] Ankyrin repeat protein, partial [Haliotis diversicolor] Thioester–containing protein–B [Azumapecten farreri] Kunitz–type proteinase inhibitor 5 II [Crassostrea gigas] Small heat shock protein 26 [Haliotis discus hannai] IgGFc–binding protein [Crassostrea gigas] Mnk [Aplysia californica] Heat shock inducible protein 70 [Haliotis diversicolor] Baculoviral IAP repeat–containing protein 7–B [Crassostrea gigas] Putative inhibitor of apoptosis [Crassostrea gigas] Baculoviral IAP repeat–containing protein 7–B [Crassostrea gigas] Nacrein [Turbo marmoratus] 60S acidic ribosomal protein P1 [Crassostrea gigas] Tripartite motif–containing protein 2 [Crassostrea gigas] Nacrein [Turbo marmoratus] Invertebrate–type lysozyme [Haliotis discus discus] Trafficking protein particle complex subunit 1 [Crassostrea gigas] Beta–ureidopropionase [Crassostrea gigas] Calmodulin [Crassostrea gigas] Ankyrin repeat protein, partial [Haliotis diversicolor] Nephrin, partial [Biomphalaria glabrata] IkB–like protein [Haliotis discus discus] Baculoviral IAP repeat–containing protein 7–B [Crassostrea gigas] GTP–binding protein Di–Ras2 [Crassostrea gigas] Tripartite motif–containing protein 45 [Crassostrea gigas] Tyrosinase–like protein tyr–2 [Pinctada martensii] Arginine kinase [Haliotis madaka] Heat shock inducible protein 70 [Haliotis diversicolor] Interferon–induced 44–like protein [Haliotis discus discus] Perlucin C protein precursor [Haliotis laevigata] Heat shock inducible protein 70 [Haliotis diversicolor] Small heat shock protein 26 [Haliotis discus hannai]
Table 1 (continued) Feature ID Guanaqueros Contig0147908 Contig0017722 Contig0032910 Contig0037229 Contig0090540 Contig0030205 Contig0127687 Contig0207515 Contig0029365 Contig0058718 Contig0051832 Contig0000560 Contig0014040 Contig0002192 Contig0023446 Contig0003505 Contig0162781 Contig0024938 Contig0031819 Contig0000811 Contig0003216 Contig0093824 Contig0011373 Contig0056561 Contig0000372 Contig0003405 Contig0002989 Contig0003068 Pta. Talca Contig0040521 Contig0018468 Contig0004869 Contig0003475 Contig0016515 Contig0084135 Contig0030041 Contig0013447 Contig0012861 Contig0004766 Contig0009998 Contig0000948 Contig0035323 Contig0008623 Contig0002667 Contig0016158 Contig0002073 Contig0005984 Contig0014022 Contig0085038 Contig0080549 Contig0022015 Contig0008481 Contig0053154 Contig0006086 Contig0002001 Contig0006069 Contig0000372 Contig0008879
Fold–change
134,1 101,4 77,8 67,2 60,5 54,1 46,4 46,1 44,4 44,0 37,2 36,3 34,6 32,1 31,6
Description
–21,5 –21,5 –24,5 –26,1 –30,7 –39,6 –73,1 –109,4 –124,2 –159,7 –189,2 –310,6 –501,3
Maspardin [Crassostrea gigas] Ribosomal protein L7, partial [Doryteuthis pealeii] Tyrosine––tRNA ligase, mitochondrial–like [Aplysia californica] Aplycalcin [Aplysia californica] Tropomyosin [Turbo cornutus] Serine/threonine–protein phosphatase 4 regulatory subunit 1 [Crassostrea gigas] Putative actin [Arion lusitanicus] MitochondrialATP synthase subunit 9 precursor–like protein [Haliotis diversicolor] Ras–related C3 botulinum toxin substrate 1 [Crassostrea gigas] Toll–like receptor 13–like [Aplysia californica] Hemocyanin isoform 1, partial [Haliotis diversicolor] Nacrein [Turbo marmoratus] Hemicentin–1 [Crassostrea gigas] Carbonic anhydrase [Haliotis tuberculata] Ferritin [Haliotis diversicolor] Heat shock protein 70 B2 [Crassostrea gigas] 40S ribosomal protein S14 [Crassostrea gigas] Ankyrin repeat domain–containing protein 50–like, partial [Aplysia californica] Temptin [Haliotis discus discus] Cysteine–rich motor neuron 1 protein–like [Aplysia californica] Perlucin A2 protein precursor [Haliotis laevigata] IkB–like protein [Haliotis discus discus] Leukocyte cell derived chemotaxin 1–like protein [Haliotis discus discus] Small heat shock protein 26 [Haliotis discus hannai] Lustrin A [Haliotis rufescens] Heat shock inducible protein 70 [Haliotis diversicolor] Putative C1q domain containing protein MgC1q23 [Mytilus galloprovincialis] 40S ribosomal protein S3a–like [Aplysia californica]
–24,9 –26,0 –26,8 –29,7 –101,9 –132,4 –234,3 –289,9 –393,8 –541,0
Carbonic anhydrase [Haliotis gigantea] COP9 signalosome complex subunit 3 [Crassostrea gigas] Neuroglian [Crassostrea gigas] Sodium– and chloride–dependent GABA transporter 2 [Crassostrea gigas] Serine/threonine–protein kinaseWNK1 [Crassostrea gigas] Vasa–like protein [Haliotis asinina] RAS and EF–hand domain–containing–like protein [Crassostrea gigas] Fibrinogen–like protein A [Crassostrea gigas] Multiple epidermal growth factor–like domains 11 [Crassostrea gigas] SCO–spondin [Crassostrea gigas] Nacre serine protease inhibitor 3 [Pinctada margaritifera] Secreted protein [Haliotis discus discus] Ankyrin [Crassostrea gigas] Angiogenic factor with G patch and FHA domains 1 [Crassostrea gigas] Low–density lipoprotein receptor–related protein 6 [Crassostrea gigas] Fibroblast growth factor receptor [Crassostrea gigas] Mu class glutathione–s–transferase [Haliotis discus discus] Sodium–dependent proline transporter [Crassostrea gigas] Fibropellin–1 [Crassostrea gigas] Vitelline coat protein 42 [Omphalius pfeifferi] Histone 2A [Aplysia californica] Small heat shock protein [Tegillarca granosa] Mesenchyme–specific cell surface glycoprotein [Crassostrea gigas] GTPase IMAP family member 1 [Crassostrea gigas] Baculoviral IAP repeat–containing protein 7–B [Crassostrea gigas] Kunitz–type proteinase inhibitor 5 II [Crassostrea gigas] Putative inhibitor of apoptosis [Crassostrea gigas] Lustrin A [Haliotis rufescens] Heat shock inducible protein 70 [Haliotis diversicolor]
201,2 75,7 67,0 60,5 59,8 56,3 53,6 53,6 50,1 47,7 46,2 46,2 45,4 45,4 44,9 42,8 42,1 41,5 40,1
transcript per million mapped reads (RPKM). The parameters considered included a minimum read length fraction = 0.8, minimum read similarity fraction = 0.9, and unspecific read match limit = 10 in relation to the reference values. Transcriptional differences were tested through clustering analysis among all T. atra populations, and between summer and winter. The metric distance was calculated using the Manhattan method, where the mean expression level in 5–6 rounds of k-means clustering was subtracted. Finally, a Kal's statistical analysis test was used to compare gene expression levels among datasets in terms of the log2 fold change (P = 0.0005; FDR corrected). During the summer, three clusters of transcription expression were mainly identified. Cluster 1 was mostly associated with Temblador, Guanaqueros y Pta. Talca populations, displaying Huasco as the population with the lowest values of gene expression. Cluster 2 evidenced an opposite pattern with a highest transcription expression from snails sampled in Huasco during the summer. An additional cluster was highly expressed in individuals from Pta. Talca population (Fig. 1AS). Gene annotation (see Table 1 and S2) revealed that most of the transcripts identified were antioxidant system and molecular chaperons such as heat shock proteins (HSPs). Herein, to evidence the transcription expression, cluster analysis of genes annotated as antioxidant and heat shock proteins from T. atra populations was performed (Fig. 1). Putative glutathione S-transferases and peroxidases evidenced high transcription
Please cite this article as: Tapia, F.J., Gallardo-Escárate, C., Spatio-temporal transcriptome analysis in the marine snail Tegula atra along centralnorthern Chile (28–31°S), Mar. Genomics (2015), http://dx.doi.org/10.1016/j.margen.2015.05.005
F.J. Tapia, C. Gallardo-Escárate / Marine Genomics xxx (2015) xxx–xxx
3
A
B
C
D
Fig. 1. Transcriptome analysis of antioxidant system and heat shock proteins from T. atra populations during the summer and winter. A, B, C and D represents the hierarchical clustering of differentially transcribed genes among Huasco, Temblador, Guanaqueros and Pta. Talca populations.
expression during the summer in Huasco-Temblador and GuanaquerosPta. Talca populations, respectively. As expected, chaperon genes such as HSP70 and 90 were highly expressed during the summer in northern populations. Interestingly, isoforms of HSPs were also expressed during the winter in all populations. As for the transcriptional patterns observed in T. atra populations during the winter, two main clusters were identified. Cluster 1 evidenced high gene expression in the Huasco and Guanaqueros populations; meanwhile cluster 2 was characterized by low transcription activity in Huasco
and Pta. Talca (Fig. 1BS). Gene annotation evidenced a wide array of expressed genes. For instance, several transcripts related to immune response and metabolisms were identified (Tables 1 and S3). Furthermore, carbonic anhydrase (CA) was also highly expressed in populations localized in Huasco and Temblador. Herein, CA is not only central to the transport and excretion of CO2 in animals (or uptake in autotrophic organisms), but is also indirectly involved in response to conditions in the ambient environment, and may reflect the metabolic rate of the organism (Zhang et al., 2012).
Please cite this article as: Tapia, F.J., Gallardo-Escárate, C., Spatio-temporal transcriptome analysis in the marine snail Tegula atra along centralnorthern Chile (28–31°S), Mar. Genomics (2015), http://dx.doi.org/10.1016/j.margen.2015.05.005
4
F.J. Tapia, C. Gallardo-Escárate / Marine Genomics xxx (2015) xxx–xxx
A
Chile
B Huasco
Temblador
Upregulated Guanaqueros Downregulated
Pta. Talca
0
500
1000
1500
2000
2500
Number of transcripts
C
Fig. 2. Climatology and differentially expressed transcripts from studied populations. A) Mean Sea Surface Temperature (SST) over the study region from MODIS-Aqua imagery for July 2009–December 2014 (left), and mean seasonal cycle of SST at coastal pixels (right). B) Number of transcripts up and downregulated (summer/winter) for each T. atra population. C) Analysis of means for fold-changes values between summer and winter for all snail populations analyzed. Shaded lines and red dots represent significant average values (α = 0.1). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)
2.4. Thermal variability and transcriptional modulation To assess the effect of thermal variability as a key factor of environmental stress, data on Sea Surface Temperature (SST) over the
study region were gathered from Level-3 images produced by the MODIS (Moderate Resolution Imaging Spectroradiometer) sensor onboard NASA's Aqua satellite. Daily images for the period July 2009– December 2014, and with a 4 × 4 km spatial resolution (N = 390)
Please cite this article as: Tapia, F.J., Gallardo-Escárate, C., Spatio-temporal transcriptome analysis in the marine snail Tegula atra along centralnorthern Chile (28–31°S), Mar. Genomics (2015), http://dx.doi.org/10.1016/j.margen.2015.05.005
F.J. Tapia, C. Gallardo-Escárate / Marine Genomics xxx (2015) xxx–xxx
were used to characterize the mean SST field and the mean annual cycle for coastal pixels along this region Fig. 2A. The coastal SST climatologies showed that Huasco, Temblador, and Guanaqueros are the warmest localities during the summer (January–March), whereas Pta. Talca is cooler both in the summer and winter months Fig. 2A. A closer examination of the SST climatology revealed that Huasco and Guanaqueros differ substantially in terms of thermal variability. Interestingly, the number of transcripts up/down-regulated was consistently lower in these populations relative to the snails sampled from Temblador and Pta. Talca (Fig. 2B). Similarly, fold changes values of up/downregulated genes showed that the individuals from Huasco and Guanaqueros display high transcriptional activities (Fig. 2C). This fact could suggest that the expression of these genes has evolved a degree of “frontloading” that potentially pre-adapts the populations to frequent heat stress and contributes to their higher thermal tolerance (Barshis et al., 2013; Gleason and Burton, 2015). These results may suggest that both upregulation and frontloading of genes may be employed as a molecular mechanism for local adaptation along the Chilean coast.
Acknowledgments This study was funded by CONICYT-Chile through FONDECYT grant 1120896 and FONDAP project 15110027.
5
Appendix A. Supplementary data Supplementary data to this article can be found online at http://dx. doi.org/10.1016/j.margen.2015.05.005. References Barshis, D., Ladner, J., Oliver, T., Seneca, F., Traylor-Knowles, N., Palumbi, S.R., 2013. Genomic basis for coral resilience to climate change. PNAS 110, 1387–1392. Conesa, A., Götz, S., García-Gómez, J.M., Terol, J., Talón, M., Robles, M., 2005. Blast2GO: a universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics 21, 3674–3676. Gleason, L.U., Burton, R.S., 2015. RNA-seq reveals regional differences in transcriptome response to heat stress in the marine snail Chlorostoma funebralis. Mol. Ecol. 24, 610–627. Haye, P.A., Segovia, N.I., Muñoz-Herrera, N., 2014. Phylogeographic structure in benthic marine invertebrates of the southeast Pacific coast of Chile with differing dispersal potential. Plos One 9, e88613. Kulikova, V., Omel'yanenko, V., 2000. Reproduction and larval development of the gastropod mollusk Tegula rustica in Peter the Great Bay, Sea of Japan. Russ. J. Mar. Biol. 26, 128–130. Nunez-Acuna, G., Tapia, F.J., Haye, P.A., Gallardo-Escarate, C., 2012. Gene expression analysis in Mytilus chilensis populations reveals local patterns associated with ocean environmental conditions. J. Exp. Mar. Biol. Ecol. 420, 56–64. Sagarin, R.D., Gaines, S.D., 2002. Geographical abundance distributions of coastal invertebrates: using one-dimensional ranges to test biogeographic hypotheses. J. Biogeogr. 29, 985–997. Valdovinos, C., 1999. Biodiversidad de moluscos chilenos: Base de datos taxonómica y distribucional. Gayana 62, 111–164. Zhang, G., Fang, X., Guo, X., 2012. The oyster genome reveals stress adaptation and complexity of shell formation. Nature 490, 49–54.
Please cite this article as: Tapia, F.J., Gallardo-Escárate, C., Spatio-temporal transcriptome analysis in the marine snail Tegula atra along centralnorthern Chile (28–31°S), Mar. Genomics (2015), http://dx.doi.org/10.1016/j.margen.2015.05.005
