Transcriptional Responses to Long-term Salinity Stress and Acclimation in the Black Mussel, Mytilus Galloprovincialis

doi:10.21203/rs.3.rs-1536907/v1

Download PDF

Article

Transcriptional Responses to Long-term Salinity Stress and Acclimation in the Black Mussel, Mytilus Galloprovincialis

https://doi.org/10.21203/rs.3.rs-1536907/v1

This work is licensed under a CC BY 4.0 License

Version 1

posted

You are reading this latest preprint version

The black mussel Mytilus galloprovincialis is an ideal organism for studying adaptation as it inhabits highly variable environments. We examined its acclimation to long-term, time-series osmostress utilizing RNA-seq based transcriptomics approaches. In our results, osmotic homeostasis was enhanced by alterations of membrane permeability and nitrogen metabolism. Cholinergic ciliary stimulation and calcium signaling were involved in the response. Genes encoding protein turnover, carbohydrate metabolism/catabolism, nucleotide metabolism, arachidonic acid metabolism, aerobic and anaerobic energy metabolisms, and apoptosis were all inversely regulated with salinity. In addition to the existing literature, we encountered two important metabolic regulations; first, involvement of anaerobic to aerobic metabolism that was suggested in bivalves before, and second, regulation of PEPCK with OXPHOS and glycolysis genes. These two regulations attracted our attention, especially in terms of their similarities to the metabolic regulations seen in cancer cells. Although there are many types/causes of cancer, there are common adaptations that support survival/proliferation, such as the aforementioned metabolic regulations and apoptosis suppression. In our study, parallel with these metabolic regulations, p53 and p63 originated apoptosis were present with the participation of TP53 apoptosis effector (PERP), TLRs, and TNFSF14. Understanding genes/pathways/modifications that trigger the mechanism that sets off cell death in this context in our study can be promising for cancer and in developing therapeutic and protective products such as vaccines.

Understanding molecular mechanisms involved in stress response is important and promoted by scientists, as it improves our knowledge about the physiology of the cell and contributes to a broad range of biological applications such as developing stress-tolerant crops or understanding the resistance of cancer cells to drugs¹. As human impact on ecosystems and climate increases, it is also essential to understand these mechanisms in predicting the course of the species in terms of their survival vs. extinction. In marine environments, it important in predicting distribution patterns, abundance, competition capacity, and ultimate fate of the species²^,³.

Sessile marine organisms that inhabit estuaries, lagoons, or intertidal and subtidal zones experience drastic fluctuations in environmental conditions as they cannot migrate to compensate for the change. They are constantly exposed to high levels of abiotic stress such as oxygen, temperature, pH and salinity change due to tides, heavy rains, groundwater discharges and other factors⁴^,³. One sessile resident of these zones is the black mussel, Mytilus galloprovincialis (Lamarck, 1819), a filter-feeder exposed to high levels of pathogens and pollution besides other varying conditions in the adult form⁵. Interestingly, although M. galloprovincialis is exposed to high levels of pathogens, it is much less susceptible to infectious outbreaks and mass mortality compared to other bivalves like oysters and clams⁵. The ability of M. galloprovincialis to adapt and survive to ever-changing local conditions and deal with environmental stress, makes it an ideal model organism for studying animal adaptation.

Salinity is one of the most important factors that affect M. galloprovincialis in many aspects such as physiology, cellular homeostasis, development, abundance, feeding, reproduction, growth, respiration, and interspecific interactions among others³^,⁶^,². M. galloprovincialis adapt to changing salinities with behavioral and molecular responses. All of the processes occur by the expression of genes that are triggered through signaling pathways to ensure organismal homeostasis⁷. While some of these changes in gene expression occur in response to salinity stress only, others respond with changes to any kind of stress an organism might encounter. In addition to osmotic alterations, salinity change can alter organisms’ metabolism, energy reserves, oxidative status, anabolic and catabolic processes, defense mechanism, DNA damage and repair response (cellular damage), apoptosis, among others, in ways related to the degree and time of exposure to stress³^,⁴^,⁶.

The effect of gradual or minor changes in environmental conditions on organisms in their habitats can be understood by long-term acclimation stress experiments⁸. This study is focused on understanding the mechanisms of long-term acclimation response of M. galloprovincialis to gradual chronic salinity change while accounting for the lower and upper limits of salinity encountered by this species in the Bosphorous and adjacent seas. In the form of time-series experiment of salinity increase, acclimation and decrease, we tried to understand the plasticity of the expression and repeatability of the response. The salinity responses of bivalves have been examined in many studies by targeting certain markers and physiological responses⁶^,⁹^,¹⁰^,¹¹^,¹²^,¹³. By examining the entire mRNA profile with RNA-seq, we tried to fill the gaps that may occur by looking at specific markers and processes. Differential gene expression (DGE) analysis, gene set enrichment analysis (GSEA), and pathway enrichment analysis (PEA) were applied to identify genes, gene ontology (GO) groups, and metabolic pathways involved in the acclimation response.

In brief, we collected mussels and applied gradual long-term salinity increase, acclimation, and decrease in a time-series experiment. Gene expression, pathway, and gene set enrichment analyzes were used with RNA-seq with de novo transcriptomics. A schematic representation of the experimental setup is given in Fig. 1.

Animals. We collected 120 adult mussels from Rumeli Hisarüstü in Istanbul Turkey (41° 5'16.60"N 29° 3'24.95"E) from a salinity of 18.7‰ and temperature of 18.1 ºC, measured using a WTW Multi 340i system. Samples were transported to the laboratory in seawater taken from the sea at the time of the sampling.

Experimental design. In the laboratory, four tanks with size 31x14x55h cm were prepared with artificial seawater (Instant Ocean). Two of the tanks were placed inside a bigger tank, Tank A, to share the same aeration pump and the cooling unit. The remaining two tanks were placed inside another big tank, Tank B, sharing the same aeration pump and the cooling unit. Bigger tanks, Tank A and Tank B, were arranged as technical replicates of each other having different aeration pumps and different cooling units. Experimental units in Tank A and Tank B were arranged as Rumeli Hisarı Treatment (RT), Rumeli Hisarı Control (RC), respectively (Fig. 1a).

Tanks were kept at 18.0 ºC temperature and baseline salinities of 18‰. Each tank had 30 mussels with approximately 800 ml of water per mussel. Each tank was constantly aerated. Mussels were fed daily (20 µl/mussel) with phytoplankton of PhytoGreen-M (Brightwell Aquatics). Every week, approximately 1/3 of the water was changed with artificial seawater having the same salinity and temperature conditions. We analyzed 48 samples with RNA-seq by pooling three samples per experimental unit for each timepoint, see below for details. The reamining samples were collected considering the chances of mortality and/or experimental problems, as backup. The amount of food and water available per mussel was kept constant by adjusting their volumes according to the number of live samples that remained in the tank, after the death of individual mussels during the experiment, when it occured.

Experimental period. The experiment lasted 84 days in total. Four timepoints, A, B, C and D, with time intervals of 21 days were determined for sampling to observe treatment effects (Fig. 1a, Fig. 1b, Fig. 1c). Control groups were kept at constant salinity of 18‰ (salinity of Bosphorus) to be used as time controls throughout the study.

The first sampling was at timepoint A, in which mussels were acclimated to baseline salinity of 18‰ for 21 days to eliminate the effect of other stressors. Next, salinity was increased by 1‰ per day until reaching 39‰ (salinity of Aegean Sea) within 21 days for the treatment groups. At timepoint B, the second sampling was conducted at 39‰ for the treatments and 18‰ for the control groups. The treatment groups were kept 21 days at 39‰ salinity for acclimation while the control groups were still at constant salinity of 18‰. After this acclimation period, sampling was conducted at timepoint C for both groups. Lastly, the salinity of the treatment groups was decreased by 1‰ per day until reaching 18‰ (within 21 days) while control groups were at constant salinity of 18‰. The last sampling at timepoint D was conducted at the end of this period (Fig. 1a, Fig. 1b).

RNA isolation and Illumina sequencing. The experiment was conducted between November 2016 and February 2017. For the sampling, gill tissue was dissected, finely cut into pieces and frozen at -80°C in DNA/RNA/Protein isolation reagent (Tri Reagent). Total RNA was extracted with spin column purification using Direct-Zol RNA MiniPrep Plus Kit following the manufacturer's protocol in March 2019. Extracted total RNAs were sent to Macrogen Inc. on dry ice for library preparation and sequencing. Quantification and quality control of total RNA was measured with 2200 TapeStation by Macrogen Inc. RNAs of three individuals for each library were pooled in equal quantities generating 16 RNA-seq libraries (Fig. 1a, Fig. 1b) to increase library complexity and decrease inter-individual variability. Eight treatment points, each with a replicate (from replicated tanks, Fig. 1a, Fig. 1b) was chosen for RNA-seq. As replicates were coming from different tanks, the tank effect was incorporated into the model matrix (Eq. 1, Eq. 2) when building the model. Library preparation was done non-stranded using TruSeq RNA Sample Prep Kit v2, following TruSeq RNA Sample Preparation v2 Guide Part # 15026495 Rev. F by Macrogen Inc. Libraries were sequenced on Illumina HiSeq2500 platform yielding 100 bp paired-end reads, ~1Gb/sample.

Data analysis. A flowchart of data analysis is depicted in Fig. 2. The processes requiring high computing power were performed at TUBITAK ULAKBIM, High Performance and Grid Computing Center (TRUBA)¹⁴. For the preparation of the raw reads, sequencing primers (SP1, SP2) and MID sequences (i5, i7) were removed with Cutadapt v1.16¹⁵ and low quality reads (Phred score < 15) and short reads (reads < 25 bp) were removed with Trimmomatic v0-2.38¹⁶. Quality was analyzed with FastQC v0.11.7¹⁷ before and after removal of adaptors, and low quality and short reads. To check for possible rRNA contamination, the amount of ribosomal RNA was examined using SortMeRNA v2.1b¹⁸.

After filtering, for comprehensive reference assembly, reads of all samples were merged. Trinity v2.8.3¹⁹ was used for constructing de novo transcriptome assembly with a default k-mer size of 25. Completeness of the transcriptome was assessed by BUSCO v3²⁰. To remove redundancy, assembled transcripts were clustered with CD-HIT-EST v4.6²¹.

Clustered assembly was queried to the following databases with an e-value of 1x10^− 7 as cut-off: blastx against NCBI NR, blastx against NCBI RefSeq protein invertebrate, blastx against Swiss-Prot, blastn against NCBI NT, blastn against NCBI RefSeq RNA, blastx against TrEMBL, blastx against NCBI RefSeq protein and blastn against RefSeq genomic. At each step results without hits were selected and queried against the next database. Blastx searches were used to identify protein products in NR, RefSeq or Swiss-Prot databases. If no protein products were found, blastn searches were conducted. Open reading frames and coding sequences were checked.

After the blast step, OmicsBox v1.3 was used for mapping and annotation²². Transcripts were assigned enzyme commission (EC) numbers by running Kyoto Encyclopedia of Genes and Genomes (KEGG)²³ services from OmicsBox²⁴. InterProScan was run through OmicsBox for annotation. Orthology assignment by eggNOG mapper was also conducted²⁵. After annotation, GSEA was done through OmicsBox²⁴. FDR was used as a rank metric in GSEA analysis. REVIGO was used for GO visualization²⁶. Summaries of enriched GO terms having FDR < 0.25 and top 50 enriched GOs by normalized enrichment score (NES) scores were interpreted.

For DGE analysis, count files were generated following the recommended workflow of StringTie²⁷ with a modification. Sample reads were mapped to de novo assembled transcriptome by HISAT2 v2.1.0²⁸. Reference annotation information was not used in the mapping process. Mapped reads were assembled using StringTie v1.3.6 for each sample. StringTie was run with --merge option for merging assembled transcripts of samples and removing redundant transcripts. Count files were generated by HTSeq v0.11.2²⁹ with sorted mapped files of each sample and merged transcripts generated in the previous step.

BAM file statistics were calculated using Samtools v1.10³⁰ -depth and -flagstat commands. The breadth of coverage was calculated as the total number of covered bases divided by the reference assembly length. The depth was calculated as the total number of bases of all mapped reads divided by reference assembly length. The mean read depth was calculated as the total number of bases of all mapped reads divided by the total number of covered bases.

A negative binomial generalized linear model was used for DGE analysis using Bioconductor package edgeR v3.26.8³¹. Counts per million (CPM) and log CPM transformations for quantifying library size differences and trimmed mean of M-values (TMM) for normalizing distributions across samples were used. Transcripts that were lowly expressed or with zero counts were removed from the transcript lists for all samples. As the differential expression condition tested for was ‘salinity stress’ between timepoints, we used ‘timepoint’ grouping factor for our experimental design, that is for a contig to be considered expressed, it should have a worthwhile expression in at least four libraries.

Contrasts examined are given in Fig. 1c. Contrast, RTvsRC was examined for examining differentially expressed (D.E.) genes between treatment and control groups. The main contrasts, RTvsRC.BA, RTvsRC.CB and RTvsRC.DC were examined for D.E. genes involved in adaptation to salinity increase from home condition, acclimation to high salinity stress, and returning to home condition, respectively. Genes that were D.E. in these main contrasts were also checked in significant lists of other contrasts, namely RTvsRC.DA, RTvsRC.DB, RTvsRC.CA and RTvsRC (Fig. 1c). In this perspective two design matrixes were constructed:

model.matrix(~ tank + treatment + timepoint + treatment:timepoint ) (1)

model.matrix(~ 0 + group + tank) (2)

model.matrix in Eq. (1) was mainly used for the comparison of the treatment group to the control group (contrast RTvsRC in Fig. 1c). Tank B, control groups and timepoint A were used as reference levels and adjustments were made accordingly. model.matrix in Eq. (2) was constructed for pairwise comparisons given in Fig. 1c. Treatment and timepoint pairs are represented as ‘group’ in the design matrix. Reference levels were the same as in Eq. (1) and adjustments were manually calculated as given in Fig. 1c.

Likelihood ratio test with tagwise dispersion was applied, as it is recommended for multi-factor experimental designs³². False discovery rate (FDR) was adjusted by Benjamini-Hochberg method³³. DGE lists were constructed with cut-off values of FDR < 0.05 and |LogFoldChange| > 2 for contrasts (Fig. 1c). KofamKOALA³⁴, EggNOG²⁵ and KAAS³⁵ mappings were added to these DGE lists for KEGG Orthology (KO).

For PEA, ComPath³⁶, which works with WikiPathways³⁷, Reactome³⁸ and KEGG³⁹ databases, were used. D.E. genes of the main contrasts (Fig. 1c.) having FDR < 0.05 were selected as genesets for the program. Enriched pathways that have adjusted p-values lower than 0.05 were sorted in descending order according to the number of genes mapped to each pathway. As M. galloprovincialis is a non-model species, orthologs genes and pathways of other species were not removed from the results except for human-specific diseases, as one of the main objectives of the study is to determine target genes and pathways and check if results support previous findings.

Limitations of the study. Due to working with a non-model organism and in order to find enriched pathways between contrasts, PEA was conducted using a set of gene symbols (significant contigs in DGE analysis) against public databases (WikiPathways, Reactome and KEGG), which must be considered with its limitations. Bulk selection of gene symbols can result in false negatives when queried against databases as different symbols can be used between databases and contig lists. Also, the results are closely linked to the data scope of each database. However, the analysis was still useful for underlining the prominent pathways, albeit with these limitations. Similarly, de novo assembly construction, mapping, and annotation for non-model organisms are also linked to data entries in the databases, which should be considered in the interpretation of the results and the context of the species.

In our study, due to budget constraints, the tanks that were designed as technical replicates for batch effects (Tank A and B) were later used as biological replicates of experimental conditions. So each replicate of an experiment condition is prepared by pooling the samples of the related tank for that condition. Since it is not possible to realize artifacts that can occur in a single tank without replicate tanks, the design was adapted in this way and the tank effect was incorporated into the model (Eq. 1, Eq. 2).

Summary of the sequencing and assembly statistics are given in Table 1 (Detailed coverage and sequencing statistics are given in Supplementary Table S1 and Table S2). 99.4% of the reads passed the quality standards for further analysis. Completeness of assembly, assessed with BUSCO analysis resulted in a score of 99.8%. Among 493,815 contigs of the assembly, 144,459 had blast hits. When eggNOG mapper results were merged, the number of GOs assigned increased from 85,837 to 480,376, and the number of enzyme codes assigned increased from 8,789 to 57,565.

Table 1

Summary of the sequencing and assembly statistics.
RNA, Sequencing and Data Quality Statistics
Average RIN value of RNAs: 9.3 (SD = 0.6) Total raw reads generated: 687,060,562 Average GC content of 16 libraries: 38.81% (SD = 0.36) Average Q20 of 16 libraries: 98.87% (SD = 0.05) Average Q30 of 16 libraries: 96.14% (SD = 0.11) Total reads left after filtering: 683,010,660
Assembly Statistics
Trinity Statistics Total trinity ‘genes’: 438,285 Total trinity transcripts: 834,295 Percent GC: 34.15 Stats based on ALL transcript contigs: Contig N10: 4712 Contig N20: 3147 Contig N30: 2277 Contig N40: 1679 Contig N50: 1213 Median contig length: 400 Average contig: 747.53 Total assembled bases: 623,656,913 Stats based on only longest isoform per ‘gene’:Contig N10: 3703 Contig N20: 2271 Contig N30: 1496 Contig N40: 1006 Contig N50: 703 Median contig length: 338 Average contig: 562.46 Total assembled bases: 246,519,298	BUSCO Result C:99.8%[S:25.1%,D:74.7%],F:0.2%,M:0.0%,n:978 976 Complete BUSCOs © 245 Complete and single-copy BUSCOs (S) 731 Complete and duplicated BUSCOs (D) 2 Fragmented BUSCOs (F) 0 Missing BUSCOs (M) 978 Total BUSCO groups searched
	CD-HIT-EST statistics Total contigs: 493,815 Average contig length: 629.40 Total assembled bases: 310,805,708 Contigs over 5 Kb in length: 2,883 Longest contig length: 36,320
Mapping and Alignment Statistics
Overall alignment rate of 16 libraries to de novo assembly: > 85% Average alignment of 16 libraries to the de novo assembly: 87.12% (SD = 0.91) Average breadth of coverage of 16 libraries: 41% (SD = 3) Average depth of 16 libraries: ~9X (SD = 0.9) Average mean read depth of 16 libraries: 21.22 (SD = 1.7)
Annotation Statistics
Contigs with blastx hit to NR: 116,225 Contigs with blastx hit to RefSeq Invertebrate: 3,485 Contigs with blastx hit to Swiss-Prot: 103 Contigs with blastn hit to NT: 20,792 Contigs with blastn hit to RefSeq RNA: 39 Contigs with blastx hit to TrEMBL: 584 Contigs with blastx hit to RefSeq Protein: 22 Contigs with blastn hit to RefSeq Genomic: 3,209 Total of blast hits: 144,459 Contigs with Blast Hits in OmicsBox: 143,999 Contigs with InterProScan: 108,887 Contigs with EggNOG: 87,900 Contigs with Mapping: 80,507 Contigs with GO Annotation: 69,399	Species#BLAST Top-Hits (Top 5) Mytilus galloprovincialis: 38,267 Mizuhopecten yessoensis:18,430 Crassostrea virginica: 13,452 Crassostrea gigas: 12,713 Bathymodiolus brooksi thiotrophic gill symbiont:1,801 Species#BLAST Hits (Top 5) Crassostrea virginica: 161,067 Mizuhopecten yessoensis: 114,775 Crassostrea gigas: 103,995 Mytilus galloprovincialis: 97,858 Pomacea canaliculata: 23,537

Gene set enrichment analysis. 417, 257, and 388 GOs were enriched (with FDR < 0.25), for salinity increase, high salinity acclimation, and salinity decrease, respectively. GO summarizations and top 50 enriched GOs by NES scores are given in Supplementary Figure S6-S17 and Supplementary Table S4. Broadly, gene enrichments were grouped under ciliary activity (enriched in salinity increase and high salinity acclimation), energy metabolism - oxidative phosphorylation (OXPHOS) (enriched in salinity decrease), gene expression regulation (in all of the treatments while ‘negative regulation of gene expression’ was enriched in salinity increase), protein metabolism (in all of the treatments, especially protein catabolic processes were highly enriched for salinity increase), nucleotide metabolism (in all treatments in REVIGO summarizations and in salinity increase and decrease in top 50 enriched GO terms), and mild DNA repair (in all treatments in REVIGO summarizations and in salinity decrease in top 50 enriched GO terms). The figures and tables referred to are given in Supplementary Table S3.

Pathway enrichment analysis. Top 50 enriched pathways (Supplementary Table S5) having FDR < 0.05 were interpreted to highlight the most prominent pathways for different salinity treatments. Pathways related to immune system, metabolism, calcium signaling, the structural integrity of tissues/cells, cellular homeostasis and processes such as growth, proliferation, etc. and their respective signaling pathways were enriched.

Pathways of immune system and innate immune system were enriched in all of the contrasts. Other enrichments that were related to immune system were pathways of pathogen recognition receptors, complement system, apoptosis and their respective signaling pathways. Pathways of NOD-like receptors (NLRs), Toll-like receptors (TLRs) and C-type lectin receptors (CLRs) were enriched in all of the three treatments, whereas RIG- like receptors (RLRs) were enriched in salinity increase and high salinity acclimation. TLRs were enriched in both salinity increase and high salinity acclimation but higher number of pathways were enriched in salinity increase. Additionally, in all of the three databases, TLRs were enriched in salinity increase and high salinity acclimation, and RLRs were enriched in salinity increase. Pathways of apoptosis were enriched in all of the contrasts, with enrichment higher in salinity increase and decrease than in high salinity acclimation. Conversely, in high salinity acclimation, pathways related to survival, such as phosphatidylinositol 3-kinase (PI3K)/protein kinase B (AKT) signaling pathway⁴⁰ and pathways of growth factors (epidermal growth factor receptor (EGFR) and fibroblast growth factor receptor (FGFR)) were enriched.

Another group of pathways that were enriched was metabolism. Pathways related to protein metabolism, and post-translational modification and transportation of proteins were highly enriched in salinity increase. For lipids, pathways related to the metabolism of fatty acids were enriched with all contrasts being prominently higher in salinity decrease. Pathways of membrane lipids that are also involved in cell-cell communication, glycosphingolipid⁴¹, and sphingolipids⁴² were enriched in salinity increase and high salinity acclimation.

Pathways related to supporting the structural integrity of tissues and cells were also enriched. These were pathways of extracellular matrix (ECM), cytoskeleton and cell-matrix adhesions. These pathways can be altered due to volume changes that act as a signal that modulates metabolism, growth and development⁴³. Mechanical stimuli are converted to chemical signals for metabolic reactions via transmembrane integrins which connect ECM and cytoskeleton at focal adhesions⁴³.

Lastly, pathways related to calcium signaling, presynaptic depolarization and calcium channel opening were enriched in salinity increase, and response to elevated platelet cytosolic Ca²⁺ was enriched in salinity decrease.

Differential gene expression analysis. There were 237,591 contigs in the libraries of count files before filtering. After zero (3,797) and lowly expressed contigs were removed, 82,620 contigs were left in the libraries. In the contrasts (Fig. 1c), total of 2944 genes, 2034 of them being unique were D.E. after cut-off was applied. Of these, the main contrasts, RTvsRC.BA, RTvsRC.CB and RTvsRC.DC had 660, 490 and 516 DEGs, respectively. DEGs grouped under energy, carbohydrate, nucleotide, lipid and aminoacid metabolisms, protein turnover, antioxidant response, ion channels/transporters and apoptosis were interpreted (Table 2–5 and Supplementary Table S6).

In this study, we examined chronic long-term stress response of M. galloprovincialis. We found that the mussels adapted to changing salinities by altering energy, carbohydrate, lipid, nitrogen, and nucleotide metabolisms. Ciliary activity and apoptosis also played important roles in adaptation to stress. Results of DGE, GSEA and pathway enrichments were compatible. In addition to these groups, GO terms related to DNA repair and gene expression regulation, and pathways related to immune system were enriched in GSEA and PEA, respectively. We generated a hypothesized pathway as given in Fig. 3.

Firstly, gill cilia were probably involved in osmosensing and commanding the cellular response^44,45. Cilia act as a pivotal osmosensor and commander, sensing and transmitting extracellular signals via receptors, effector proteins and transcription factors to cells for the response⁴⁴. In the study, GO terms related with ciliary activity were enriched in salinity increase and high salinity acclimation. Contigs of cholinergic ciliary stimulation that are linked to membrane depolarization, cytosolic calcium increase and cell-to-cell signaling were D.E.^46,47,48,45. These were acetylcholine receptors^49,50,51 (Table 4), voltage-gated calcium channels⁴⁶ (Table 4), dyneins⁵² (Table 3), and a G-protein-coupled receptors (GPRC), cadherin EGF LAG seven-pass G-type receptor 3 (CELSR3)⁵³ (Table 3). Those contigs that were involved in calcium regulation and signaling were transient receptor potential melastatin channels (TRPMs)⁵⁴ (Table 4), voltage-gated calcium channels⁴⁶ (Table 4) and stanniocalcin⁵⁵ (Table 2).

Secondly, osmotic homeostatis was likely achieved. Since M. galloprovincialis is an osmoconformer, osmotic homeostasis and volume control in salinity changes depend solely on isosmotic intracellular regulation, and are mainly achieved by regulation of cell membrane permeability and intracellular osmolytes⁴. In the study, osmolyte/cell volume regulation was mainly enhanced by changing the permeability of water and ions across the cell membrane⁵⁶ and by alteration of nitrogen metabolism⁹.

Permeability and fluidity of the membrane are altered by the change in lipid composition of the membrane which triggers signaling cascades that lead to osmolyte regulation^57,58,59. In the study, change in lipid composition of the membrane by arachidonic acid metabolism and polyketide synthase probably led to osmolyte/volume regulation (Table 2). For example, regulation of sodium-potassium (Na⁺/K⁺) pump (Table 4) in the study was in line with alteration of membrane lipid composition with arachidonic acid metabolism (Table 2), as eicosanoids are modulators of this pump⁶⁰. Our findings match with the literature in that arachidonic acid metabolism is involved in a signaling cascade that leads to osmolyte release triggered by cell swelling^61,62, and Na⁺/K⁺ pump is involved in volume regulation due to the Donnan effect and swelling tendency caused by osmotic pressure⁶³. In our study, genes involved in the process were inversely regulated with salinity (downregulated with salinity increase and upregulated with decrease; cytosolic phospholipase A2 (cPLA₂), arachidonate 5-lipoxygenase (ALOX5), arachidonate 12-lipoxygenase (ALOX12B), cytochrome P450 family 2 subfamily J (CYP2J), prostaglandin-H2 D-isomerase / glutathione transferase (HPGDS), polyketide synthase (pksJ)⁵⁷ (Table 2), cholecystokinin receptor (CCKAR)^64,65,66,67 (Table 3), annexin (ANXA)^68,69 (Table 2), stanniocalcin⁵⁵ (Table 2), Na⁺/K⁺ transporting, ATPase subunit beta 1 (ATP1B)⁵⁷ (Table 4). Regulation of cytosolic Ca^+ 2 and calcium signaling in the study probably had implications for the process^62,58,70,46, as well.

Nitrogen metabolism was also altered in the study probably for osmolyte/cell volume regulation. Under hyperosmotic stress, intracellular organic osmolytes such as alanine, betaine (glycine betaine), proline, and glycine are accumulated to maintain cell volume⁷¹. In bivalves, serine, proline, and glycine are among the most common amino acids involved in osmoregulation^71,72. In this study, matching with this literature, genes involved in the synthesis of these amino acids were upregulated with salinity increase (contgis mapped to D-3-phosphoglycerate dehydrogenase (serA), phosphoserine aminotransferase (serC), glycine hydroxymethyltransferase (glyA), pyrroline-5-carboxylate reductase (proC)) (Table 2). Additionally, the regulation of glycine transporters (Table 2, Table 4) in the study was also consistent with the probable usage of glycine as an osmolyte for osmoregulation. Solute carrier family 6 member 5 (SLC6A5) transports glycine into the cytosol and was upregulated with salinity increase. This transporter has dual roles and is also involved in neuronal excitability besides osmotic stress response^73,62. The regulation in the study was in the opposite direction to that observed in Meng et al. (2013)⁶²; hypo-osmotic adaptation of Crassostrea gigas. Another glycine transporter (Solute carrier family 25 member 38 (SLC25A38)) that belongs to SLC25 family of mitochondrial transporters used in energy conversion and cell maintenance⁷⁴ was inversely regulated with salinity and in the same direction with energy demand (downregulated with salinity increase).

Additionally, contigs that are intermediates of proline synthesis, delta-1-pyrroline-5-carboxylate synthetase (P5CS) and ornithine aminotransferase (OAT), were upregulated with salinity decrease (Table 2), showing the possibility of involvement of proline in responding to chronic salinity decrease as well. A more probable scenario is that ornithine and α-ketoglutarate can form with OAT, rather than synthesis of proline with proC that can also feed the TCA cycle for energy production (α-ketoglutarate) at low salinities. Another gene, choline dehydrogenase (betA), involved in synthesis of another organic osmolyte, betaine (glycine betaine) was also upregulated with salinity decrease (Table 2). In bivalves and also in other species betaine is accumulated with hyperosmotic stress as an osmoprotectant⁷¹. Due to regulation in the study, other functions likely determined its regulation more so than the osmoprotective function such as its function with high energy metabolism (as low salinity is high energy expenditure in Mytilus), as further metabolization of choline may supply electron carriers (NADH) and ubiquinol for oxidative phosphorylation⁷⁵.

Finally, in terms of amino acid metabolism, carnosine synthase (CRNS1) was downregulated with salinity increase and upregulated with salinity decrease (Table 2) in agreement with Meng et al. (2013)⁶² which showed involvement of beta-alanine metabolism in hypoosmotic stress in the oyster C. gigas. CRNS1 regulation in our study can be related to carnosine’s buffering capacity against acidification⁷⁶ that can be caused by L-lactate dehydrogenase (LDHA) and V-type H⁺-transporting ATPase (ATPV0C) (Table 2). Its regulation was also coherent with arachidonic acid and aerobic metabolisms (Table 2) that can cause lipid peroxidation or reactive oxygen species (ROS) generation⁷⁷.

Protein turnover was also higher at lower salinity, represented by regulation of fucosidase (FUCA, involved in glycoprotein degradation⁷⁸), cathepsin B, L (participate in protein turnover⁷⁹), proteoglycans and ribosomal proteins (involved in protein synthesis and ribosomal biogenesis⁸⁰); all downregulated with salinity increase (Table 3). GOs related to protein catabolic processes were also highly enriched for salinity increase in GSEA analysis (Supplementary Fig. S7, S8, S9 and Table S4). Amino acids derived from the degradation of proteins can be recycled to synthesize new proteins, converted to other amino acids or used in energy metabolism by converting to TCA cycle intermediates or glucose. In this study, nitrogen metabolism was likely activated by cellular needs driven by osmoregulation, energy metabolism, and the requirements of the salinity conditions.

High salinity is associated with low energy expenditure and vice versa in Mytilus species^6,11,13,12. Regulation of cell volume by nitrogen metabolism in bivalves increases demands of oxidative metabolism, mobilization of reserves, and oxygen consumption rates⁹. Protein synthesis is one of the most energy-consuming processes in the cell⁸¹, and in the study, high protein turnover probably increased the energy demand of low salinity. Additionally, ionic regulation such as Na⁺/K⁺ pump that requires a large portion of the cell’s energy^56,82,83 also caused increase in energy demand (at low salinity/steady-state). Altered nitrogen metabolism and ion channel activity for osmoregulation probably reflected on energy production and metabolism. M. galloprovincialis belongs to the family Mytilidae that includes euryoxic and facultative anaerobes⁸⁴. In this study, alteration of energy metabolism was observed both in the genes of anaerobic and aerobic metabolisms (Table 2). For aerobic metabolism, oxidative phosphorylation contigs were downregulated with salinity increase (Table 2), matching with the literature as low salinity is associated with high energy expenditure in Mytilus species^6,11,13,12.

Remarkably, a contig of the anaerobic metabolism, specifically LDHA, was also downregulated with salinity increase (Table 2), in parallel with the regulation of OXPHOS contigs (anaerobic metabolism was also accompanied with monocarboxylate transporters (SLC16 in Table 4) in the study). Normally, the production of lactate from pyruvate happens when there is not sufficient oxygen⁸⁵. In this study, no indication of oxygen deprivation was found, such as differential expression of hypoxia-inducible factors (HIFs). Beyond oxygen levels, in bivalves, switching to anaerobic metabolism coupled with depressed metabolism is an energy-saving strategy in response to environmental stressors like salinity changes⁸⁴. Likewise, in our study, the involvement of LDHA could be due to salinity change. However, rather than shifting to anaerobic metabolism, genes of aerobic metabolism were also regulated parallelly (with other metabolism genes) also pointing to high metabolic rates at steady-state compared to chronic salinity increased state. Therefore metabolic depression and shift to anaerobic metabolism is not thought to be involved. Parallel regulation of aerobic (OXPHOS) and anaerobic (LDHA) metabolism contigs in the study could help in providing additional energy to the cells⁸⁶; and involvement of functional anaerobism when energy production of aerobiosis is not sufficient was suggested previously in bivalves⁸⁴. Besides, during conversion of pyruvate to lactate, NAD⁺ necessary for glycolysis is also regenerated.

This regulation was accompanied by a high glycolytic rate (glyceraldehyde-3-phosphate dehydrogenase (GAPDH) was downregulated with salinity increase also (Table 2)). High expression of GAPDH and LDHA at timepoint A, when compared to B could also potentially indicate the Pasteur effect, which refers to the inverse relation between O₂ consumption and glycolytic rate⁸⁷. However, due to two reasons, this probably was not the case: first, aerobic metabolism genes were also highly expressed, indicating the involvement of anaerobic metabolism as an additional energy supply to aerobic metabolism, as stated above. And second, bivalves can depress metabolism to decrease energy demand by decreasing energy expenditure rather than increasing glycolytic rate to compensate ATP deficiency of anaerobic metabolism, and therefore, can reverse the Pasteur effect^87,88. In the study, carbohydrate catabolism was also in agreement with the regulation of energy metabolism (FUCA, lactase-phlorizin hydrolase (LCT), Chitinase were also inversely regulated with salinity (Table 2)), as carbohydrates are the priority respiratory substrates in bivalves⁸⁹.

Mitochondrial uncoupling protein, UCP4 was parallelly regulated (Table 2, Table 4) with oxidative phosphorylation contigs (downregulated with salinity increase) which decrease ROS generation with little impact on membrane potential and ATP generation^90,91,92. To assess oxidative status, antioxidant enzymes were also examined (Table 3). Antioxidant expression was in opposite direction with salinity; Glutathione s-transferases (GSTs) contigs were downregulated with salinity increase (Table 3). In the study, phospholipase activity in arachidonic acid metabolism causes hydrolysis of phospholipids that can further result in the generation of eicosanoids, and lipid peroxidation by the enzymes lipoxygenase (LOX) and cyclooxygenase (COX)⁹³. GST downregulation with salinity increase (Table 3), parallel to both phospholipase activity (Table 2) and aerobic metabolism (Table 2), could therefore be a response to both protect from ROS due to energy metabolism, and lipid peroxidation due to COX and LOX (Table 2). DNA damage response was also examined and no significant D.E. genes were found except ATP-dependent DNA helicase PIF1 (downregulated with salinity decrease) (Table 3). Although we could not find a noteworthy DNA damage response in DGE analysis, cellular response to DNA damage stimulus (GO:0006974) and DNA repair (GO:0006281) was enriched in salinity decrease in top 50 enriched list of GSEA analysis (Supplementary Table S4).

Finally related with nucleotide metabolism, ribose-phosphate pyrophosphokinase (PRPS) and cytidine deaminase (cdd) were inversely regulated with salinity (Table 2), indicating altered nucleotide metabolism at lower salinities. PRPS downregulation indicates attenuated nucleotide synthesis with salinity increase⁹⁴ and cdd upregulation with salinity decrease indicates higher pyrimidine recycling and nucleic acid synthesis at lower salinities⁹⁵.

Besides metabolism, salinity change also alters inflammatory and immune responses, affecting lysozyme activity, humoral responses, immune-related genes, and apoptosis³. For example, alteration of membrane phospholipids (Table 2) and further metabolization of arachidonic acid leads to eicosanoids that can bind to their receptors and mediate immune and inflammatory responses. We observed regulation of genes related to apoptosis, TP53 apoptosis effector (PERP), TLRs and tumor necrosis factor ligand superfamily member 14 (TNFSF14) in our experiment (Table 5). These genes were downregulated with salinity increase, which may also be related to energy need and altered metabolism at lower salinities compared to higher salinities. Regulation of PERP also indicated p53 and p63 originated apoptosis. Inhibitor of apoptosis (IAPs) regulation probably increased stress tolerance in the study. Multiple epidermal growth factor-like domains protein 10 (MEGF10) was inversely regulated with salinity, downregulated with salinity increase and upregulated with salinity decrease (Table 5). Regulation of MEGF10 supported regulation of PERP and apoptosis, as it is involved in the clearance of apoptotic cells. It is also negative regulator of cell cycle affecting cell motility⁹⁶. PEA (Supplementary Table S5) also supported this result: pathways of apoptosis were enriched in all of the contrasts, with enrichment higher in salinity increase and decrease, rather than during the (high salinity) acclimation period. Conversely, during acclimation, pathways related to survival such as PI3K-AKT signaling pathway that is activated by cellular stimuli and regulates transcription, translation, proliferation, growth, and survival⁴⁰, and pathways of growth factors (EGFR, FGFR) were enriched. The result supported importance of cell cycle control in osmotic stress response in Mytilus mussels, as suggested by Lockwood and Somero (2011)⁹⁷.

To sum up, we detected regulation of cholinergic ciliary stimulation in osmostress adaptation. Genes of protein turnover (Table 3), carbohydrate metabolism/catabolism, nucleotide metabolism, arachidonic acid metabolism, aerobic and anaerobic energy metabolisms (Table 2) and apoptosis (Table 5) were all inversely regulated with salinity. In addition to the existing literature, we encountered two important metabolic regulations; first, the involvement of functional anaerobism that was suggested in bivalves before and second, regulation of phosphoenolpyruvate carboxykinase (PEPCK) with OXPHOS and glycolysis genes (Table 2). These two regulations attracted our attention especially in terms of their similarities to the metabolic regulations seen in cancer cells.

Cancer cells can experience hypoxia due to lack of capillary network for sufficient O₂ delivery⁹⁸. However, in cancer cells, lactate is produced (with high rate of glucose uptake and glycolysis) even in the presence of oxygen (aerobic glycolysis), the phenomenon being referred to as Warburg effect⁸⁵. In many cancer cells, respiration is retained with Warburg effect indicating glycolytic reprogramming, rather than dysfunctioning of the respiration⁹⁹. As ATP production is inefficient in fermentation compared to respiration, several explanations are made with regards to Warburg effect¹⁰⁰. One of the explanations involves feeding of the branching pathways from glycolysis to ensure macromolecular synthesis that is needed for the proliferation of tumour cells⁸⁵; in Warburg effect, metabolic pathways are enhanced, which supports macromolecular synthesis, such as pentose phosphate pathway, nucleotide biosynthesis, amino acid biosynthesis and lipid homeostasis^85,101. In the study, LDHA (anaerobic metabolism) and oxidative phosphorylation (aerobic metabolism) genes were parallelly downregulated with salinity increase (Table 2). The regulation was a Warburg effect-like regulation of vertebrates with a profile of respiratory gene expression that is also seen in some cancer tissues^99,85,102.

In our study, the second metabolic regulation, which was similar to the metabolic regulation seen in cancer cells, was the parallel regulation of genes of reciprocally regulated processes; glycolysis and gluconeogenesis. Normally these are reciprocal processes regulated by allosteric enzymes according to the energy charge of the cell. In our study, contigs of these processes, GAPDH and PEPCK were both downregulated with salinity increase (Table 2). Parallel regulation of these processes is also seen in some cancer types that rely on oxidative phosphorylation in limited glucose supply, where PEPCK is regulated to replenish glycolysis intermediates by utilizing TCA metabolites and carbon metabolites, as alternative carbon sources^103,104. Likewise, in lung cancer cells, lactate is converted to phosphoenolpyruvate with PEPCK in glucose deprivation¹⁰⁴. Regulation of PEPCK was found to be mediated by P53 via NAD⁺-dependent protein deacetylase sirtuin 6 (SIRT6) in cancer metabolism^105,104. In our study, SIRT6 was also upregulated with salinity increase (Table 2), which in return probably downregulates PEPCK (maybe there was no more need to utilize TCA metabolites) and decreases glucose uptake as high salinity is a relatively low energy expenditure state, when compared to low salinity, as mentioned above.

Cancer cells exhibit deregulated proliferation and survival mechanisms such as depressing cell death¹⁰⁶, and these regulations (Warburg effect and PEPCK) in cancer provide metabolites and ensure macromolecular synthesis that is needed for the uncontrolled proliferation of tumour cells. Likewise, Warburg effect-like regulation with respiration (contigs of oxidative phosphorylation and LDHA) and parallel regulation of PEPCK (contig of gluconeogenesis) with GAPDH in our study (Table 2) probably resulted from the need for energy and metabolites required to respond to low salinity stress (steady-state). As another important point, parallelly regulated, we observed apoptosis initiation with the participation of TP53 apoptosis effector (PERP), TLRs and TNFSF14 in our experiment (Table 5). Regulation of PERP also indicated p53 and p63 originated apoptosis. In the study, induction of apoptosis and IAPs (Table 5) probably managed and enhanced adaptation and durability to stress, by balancing between survival and death, caused by the salinity change.

The result of our study indicated M. galloprovincialis as a candidate species for studying metabolic regulation and pathways that may contribute to our understanding of cancer. Especially mechanism that sets off cell death and metabolic regulations seen (in cancer) can be related/coupled. Understanding the genes and pathways that trigger cell death under these conditions in this study and identifying RNA modifications that cause difference can be promising for cancer research, such as in developing therapeutic and protective products like vaccines.

Plasticity of the Response. In DGE analysis, contigs of oxidative phosphorylation genes were found to be downregulated with salinity increase, in line with high glycolytic rate and energy expenditure of low salinity (Table 2). The related contigs were not D.E. in salinity decrease, however overall the same trend related with high energy expenditure of low salinity was observed (FUCA and LCT were upregulated with salinity decrease) (Table 2). The reason was probably related to the acclimation period, as salinity increase was after the acclimation to steady-state condition, and mussels were more acclimated to low salinity whereas mussels were acclimated to high salinity before salinity decrease. On the other hand, in GSEA the result was vice versa; the same processes were also enriched in salinity decrease (for OXPHOS, Supplementary Fig. S15, S16, S17 and Table S4) verifying plasticity of the regulation and results of DGE. The same patterns were also seen in ribosomal proteins (Table 3 and Supplementary Table S4). This was probably due to the utilization of whole gene lists in GSEA but a subset in DGE analysis (i.e. those which have FDR < 0.05), with DGE analysis providing the most important hits.

Acknowledgements

The authors would like to thank Dr. Andrzej Furman and Dr. Ulaş Tezel for their valuable advice and help throughout the study. E.I. would also like to thank an anonymous supporter whose support was invaluable throughout the study. The financial support provided by the Boğaziçi University Research Fund grant (17Y00D1) is greatly acknowledged. The numerical calculations were partially performed at TUBITAK ULAKBIM, High Performance and Grid Computing Center (TRUBA resources) which is also greatly acknowledged.

Author contributions

E.I. Conceptualization, experimental design, data analysis, discussion and interpretation of the results, writing the paper; B.Z.H. supervision and contribution to data analysis; I.C. processing of the samples; R.B. conceptualization, supervision, funding acquisition and editing.

Data availability

RNA-seq reads of pooled libraries are available at the National Center for Biotechnology Information SRA database under the accession PRJNA824625 (https://www.ncbi.nlm.nih.gov/sra/PRJNA824625). The data that supports the findings of this study are available within the article and supplementary information. Additional data can be provided from the corresponding author upon request.

Additional Information

The authors declare no competing interests.

Grover, A. Molecular biology of stress responses. Cell Stress Chaperones 7, 1 (2002).
Fields, P. A., Zuzow, M. J. & Tomanek, L. Proteomic responses of blue mussel (Mytilus) congeners to temperature acclimation. J. Exp. Biol. 215, 1106–1116 (2012).
Pourmozaffar, S. et al. The role of salinity in physiological responses of bivalves. Rev. Aquac. 12, 1548–1566 (2020).
Rivera-Ingraham, G. A. & Lignot, J. H. Osmoregulation, bioenergetics and oxidative stress in coastal marine invertebrates: raising the questions for future research. J. Exp. Biol. 220, 1749–1760 (2017).
Figueras, A., Moreira, R., Sendra, M. & Novoa, B. Genomics and immunity of the Mediterranean mussel Mytilus galloprovincialis in a changing environment. Fish Shellfish Immunol. 90, 440–445 (2019).
Freitas, R. et al. Effects of seawater acidification and salinity alterations on metabolic, osmoregulation and oxidative stress markers in Mytilus galloprovincialis. Ecol. Indic. 79, 54–62 (2017).
De Nadal, E., Ammerer, G. & Posas, F. Controlling gene expression in response to stress. Nat. Rev. Genet. 2011 1212 12, 833–845 (2011).
Anestis, A., Lazou, A., Pörtner, H. O. & Michaelidis, B. Behavioral, metabolic, and molecular stress responses of marine bivalve Mytilus galloprovincialis during long-term acclimation at increasing ambient temperature. Am. J. Physiol. - Regul. Integr. Comp. Physiol. 293, 911–921 (2007).
Peteiro, L. G. et al. Responses to salinity stress in bivalves: Evidence of ontogenetic changes in energetic physiology on Cerastoderma edule. Sci. Reports 2018 81 8, 1–9 (2018).
Hamer, B. et al. Effect of hypoosmotic stress by low salinity acclimation of Mediterranean mussels Mytilus galloprovincialis on biological parameters used for pollution assessment. Aquat. Toxicol. 89, 137–151 (2008).
Riisgård, H. U., Lüskow, F., Pleissner, D., Lundgreen, K. & López, M. Á. P. Effect of salinity on filtration rates of mussels Mytilus edulis with special emphasis on dwarfed mussels from the low-saline Central Baltic Sea. Helgol. Mar. Res. 67, 591–598 (2013).
Sanders, T., Schmittmann, L., Nascimento-Schulze, J. C. & Melzner, F. High calcification costs limit mussel growth at low salinity. Front. Mar. Sci. 5, 352 (2018).
Tedengren, M. & Kautsky, N. Comparative stress response to diesel oil and salinity changes of the blue mussel, Mytilus edulis from the baltic and north seas. Ophelia 28, 1–9 (1987).
TRUBA. Turkish National e-Science e-Infrastructure. https://www.truba.gov.tr.
Martin, M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet.journal 17, 10–12 (2011).
Bolger, A. M., Lohse, M. & Usadel, B. Trimmomatic: A flexible trimmer for Illumina sequence data. Bioinformatics 30, 2114–2120 (2014).
Andrews, S. & others. FastQC: a quality control tool for high throughput sequence data. 2010. Https://Www.Bioinformatics.Babraham.Ac.Uk/Projects/Fastqc/ vol. 5 http://www.bioinformatics.babraham.ac.uk/projects/ http://www.bioinformatics.babraham.ac.uk/projects/fastqc/ (2010).
Kopylova, E., Noé, L. & Touzet, H. SortMeRNA: Fast and accurate filtering of ribosomal RNAs in metatranscriptomic data. Bioinformatics 28, 3211–3217 (2012).
Grabherr, M. G. et al. Trinity: reconstructing a full-length transcriptome without a genome from RNA-Seq data. Nat. Biotechnol. 29, 644 (2011).
Seppey, M., Manni, M. & Zdobnov, E. M. BUSCO: Assessing Genome Assembly and Annotation Completeness. Methods Mol. Biol. 1962, 227–245 (2019).
Li, W. & Godzik, A. Cd-hit: a fast program for clustering and comparing large sets of protein or nucleotide sequences. Bioinformatics 22, 1658–1659 (2006).
Götz, S. et al. High-throughput functional annotation and data mining with the Blast2GO suite. Nucleic Acids Res. 36, 3420–3435 (2008).
Ogata, H. et al. KEGG: Kyoto Encyclopedia of Genes and Genomes. Nucleic Acids Res. 27, 29–34 (1999).
OmicsBox – Bioinformatics Made Easy, BioBam Bioinformatics. https://www.biobam.com/omicsbox (2019).
Huerta-Cepas, J. et al. eggNOG 5.0: a hierarchical, functionally and phylogenetically annotated orthology resource based on 5090 organisms and 2502 viruses. Nucleic Acids Res. 47, D309–D314 (2019).
Supek, F., Bošnjak, M., Škunca, N. & Šmuc, T. REVIGO Summarizes and Visualizes Long Lists of Gene Ontology Terms. PLoS One 6, e21800 (2011).
Pertea, M. et al. StringTie enables improved reconstruction of a transcriptome from RNA-seq reads. Nat. Biotechnol. 2015 333 33, 290–295 (2015).
Kim, D., Langmead, B. & Salzberg, S. L. HISAT: a fast spliced aligner with low memory requirements. Nat. Methods 2015 124 12, 357–360 (2015).
Anders, S., Pyl, P. T. & Huber, W. HTSeq—a Python framework to work with high-throughput sequencing data. Bioinformatics 31, 166–169 (2015).
Li, H. et al. The Sequence Alignment/Map format and SAMtools. Bioinformatics 25, 2078–2079 (2009).
Robinson, M. D., McCarthy, D. J. & Smyth, G. K. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics 26, 139–140 (2010).
Chen, Y. et al. edgeR: differential analysis of sequence read count data User’s Guide.
Benjamini, Y. & Hochberg, Y. Controlling the False Discovery Rate: A Practical and Powerful Approach to Multiple Testing. J. R. Stat. Soc. Ser. B 57, 289–300 (1995).
Aramaki, T. et al. KofamKOALA: KEGG Ortholog assignment based on profile HMM and adaptive score threshold. Bioinformatics 36, 2251–2252 (2020).
Moriya, Y., Itoh, M., Okuda, S., Yoshizawa, A. C. & Kanehisa, M. KAAS: an automatic genome annotation and pathway reconstruction server. Nucleic Acids Res. 35, W182–W185 (2007).
Domingo-Fernández, D., Hoyt, C. T., Bobis-Álvarez, C., Marín-Llaó, J. & Hofmann-Apitius, M. ComPath: an ecosystem for exploring, analyzing, and curating mappings across pathway databases. npj Syst. Biol. Appl. 2018 41 4, 1–8 (2018).
Slenter, D. N. et al. WikiPathways: a multifaceted pathway database bridging metabolomics to other omics research. Nucleic Acids Res. 46, D661–D667 (2018).
Fabregat, A. et al. The Reactome Pathway Knowledgebase. Nucleic Acids Res. 46, D649–D655 (2018).
Kanehisa, M., Furumichi, M., Tanabe, M., Sato, Y. & Morishima, K. KEGG: new perspectives on genomes, pathways, diseases and drugs. Nucleic Acids Res. 45, D353–D361 (2017).
Janku, F., Yap, T. A. & Meric-Bernstam, F. Targeting the PI3K pathway in cancer: are we making headway? Nat. Rev. Clin. Oncol. 2018 155 15, 273–291 (2018).
D’Angelo, G., Capasso, S., Sticco, L. & Russo, D. Glycosphingolipids: synthesis and functions. FEBS J. 280, 6338–6353 (2013).
Merrill, A. H. Sphingolipids. in Biochemistry of Lipids, Lipoproteins and Membranes (ed. Vance, D. E. & Vance, C. E.) 363–397 (Elsevier, 2008).
Low, S. Y. & Taylor, P. M. Integrin and cytoskeletal involvement in signalling cell volume changes to glutamine transport in rat skeletal muscle. J. Physiol. 512, 481 (1998).
Pedersen, L. B., Schrøder, J. M., Satir, P. & Christensen, S. T. The ciliary cytoskeleton. Compr. Physiol. 2, 779–803 (2012).
Adams, M. The primary cilium: an orphan organelle finds a home. Nat. Educ. 3, 54 (2010).
Clapham, D. E. Calcium Signaling. Cell 131, 1047–1058 (2007).
Zagoory, O., Braiman, A. & Priel, Z. The mechanism of ciliary stimulation by acetylcholine: Roles of calcium, PKA, and PKG. J. Gen. Physiol. 119, 329–339 (2002).
Catapane, E. J., Nelson, M., Adams, T. & Carroll, M. A. Innervation of Gill Lateral Cells in the Bivalve Mollusc Crassostrea virginica Affects Cellular Membrane Potential and Cilia Activity. J. Pharmacol. reports 1, (2016).
BURBRING, E., BURN, J. H. & SHELLEY, H. J. Acetylcholine and ciliary movement in the gill plates of Mytilus edulis. Proc. R. Soc. London. Ser. B - Biol. Sci. 141, 445–466 (1953).
Aiello, E. L. Factors Affecting Ciliary Activity on the Gill of the Mussel Mytilus edulis. Physiol. Zool. 33, 120–135 (1960).
Aiello, E. & Paparo, A. A role for acetylcholine in the regulation of ciliary activity. Comp. Gen. Pharmacol. 5, 285–297 (1974).
King, S. M. Composition and assembly of axonemal dyneins. Dyneins Biol. Dynein Mot. Second Ed. 163–201 (2017) doi:10.1016/B978-0-12-809471-6.00005-X.
Uniprot. https://www.uniprot.org/uniprot/Q91ZI0.
Fernández, J. A. et al. Voltage-and cold-dependent gating of single TRPM8 ion channels. J. Gen. Physiol. 137, 173–195 (2011).
Roch, G. J. & Sherwood, N. M. Stanniocalcin Has Deep Evolutionary Roots in Eukaryotes. Genome Biol. Evol. 3, 284–294 (2011).
Deaton, L. Osmotic and Ionic Regulation: Cells and Animals. in Osmotic and Ionic Regulation: Cells and Animals (ed. Evans, D. H.) 107–135 (CRC Press, 2008).
Bhoite, S. & Roy, R. Role of membrane lipid in osmoregulatory processes during salinity adaptation: A study with chloride cell of mud crab, Scylla serrata. Mar. Freshw. Behav. Physiol. 46, 287–300 (2013).
Ortiz, G. G. et al. Oxidative Stress: Love and Hate History in Central Nervous System. Adv. Protein Chem. Struct. Biol. 108, 1–31 (2017).
Nelson, J. et al. Relationship between membrane permeability and specificity of human secretory phospholipase A2 isoforms during cell death. Biochim. Biophys. Acta - Biomembr. 1808, 1913–1920 (2011).
Therien, A. G. & Blostein, R. Mechanisms of sodium pump regulation. Am. J. Physiol. - Cell Physiol. 279, (2000).
Lambert, I. H. Regulation of the Cellular Content of the Organic Osmolyte Taurine in Mammalian Cells. Neurochem. Res. 2004 291 29, 27–63 (2004).
Meng, J. et al. Genome and Transcriptome Analyses Provide Insight into the Euryhaline Adaptation Mechanism of Crassostrea gigas. PLoS One 8, e58563 (2013).
Mongin, A. A. & Orlov, S. N. Mechanisms of cell volume regulation and possible nature of the cell volume sensor. Pathophysiology 8, 77–88 (2001).
Lankisch, T. O., Nozu, F., Owyang, C. & Tsunoda, Y. High-affinity cholecystokinin type A receptor/cytosolic phospholipase A2 pathways mediate Ca2 + oscillations via a positive feedback regulation by calmodulin kinase in pancreatic acini. Eur. J. Cell Biol. 78, 632–641 (1999).
González, A. et al. Cholecystokinin-Evoked Ca2 + Waves in Isolated Mouse Pancreatic Acinar Cells Are Modulated by Activation of Cytosolic Phospholipase A2, Phospholipase D, and Protein Kinase C. Biochem. Biophys. Res. Commun. 261, 726–733 (1999).
Kano, M. et al. Increased expression of gallbladder cholecystokinin: A receptor in prairie dogs fed a high-cholesterol diet and its dissociation with decreased contractility in response to cholecystokinin. J. Lab. Clin. Med. 139, 285–294 (2002).
Tsunoda, Y. & Owyang, C. High-affinity CCK receptors are coupled to phospholipase A2 pathways to mediate pancreatic amylase secretion. Am. J. Physiol. - Gastrointest. Liver Physiol. 269, (1995).
Demonbreun, A. R. et al. Recombinant annexin A6 promotes membrane repair and protects against muscle injury. J. Clin. Invest. 129, 4657–4670 (2019).
Tao, M., Isas, J. M. & Langen, R. Annexin B12 Trimer Formation is Governed by a Network of Protein-Protein and Protein-Lipid Interactions. Sci. Reports 2020 101 10, 1–11 (2020).
Evans, J. H., Gerber, S. H., Murray, D. & Leslie, C. C. The Calcium Binding Loops of the Cytosolic Phospholipase A2 C2 Domain Specify Targeting to Golgi and ER in Live Cells. Mol. Biol. Cell 15, 371–383 (2004).
Deaton, L. E. Hyperosmotic volume regulation in the gills of the ribbed mussel, Geukensia demissa: rapid accumulation of betaine and alanine. J. Exp. Mar. Bio. Ecol. 260, 185–197 (2001).
Hanson, J. A. & Dietz, T. H. The role of free amino acids in cellular osmoregulation in the freshwater bivalve Ligumia subrostrata (Say). Can. J. Zool. 54, 1927–1931 (1976).
Huang, X. et al. The Drosophila inebriated-Encoded Neurotransmitter/Osmolyte Transporter: Dual Roles in the Control of Neuronal Excitability and the Osmotic Stress Response. Genetics 160, 561–569 (2002).
Ruprecht, J. J. & Kunji, E. R. S. The SLC25 Mitochondrial Carrier Family: Structure and Mechanism. Trends Biochem. Sci. 45, 244–258 (2020).
Mailloux, R. J. et al. Choline and dimethylglycine produce superoxide/hydrogen peroxide from the electron transport chain in liver mitochondria. FEBS Lett. 590, 4318–4328 (2016).
Doğru-Abbasoğlu, S., Koçak-Toker, N. & Uysal, M. Carnosine as a Putative Antioxidant in Usage Against Liver Disease. Liver Oxidative Stress Diet. Antioxidants 295–304 (2018) doi:10.1016/B978-0-12-803951-9.00024-0.
Szabadfi, K., Pinter, E., Reglodi, D. & Gabriel, R. Neuropeptides, Trophic Factors, and Other Substances Providing Morphofunctional and Metabolic Protection in Experimental Models of Diabetic Retinopathy. Int. Rev. Cell Mol. Biol. 311, 1–121 (2014).
Artigaud, S. et al. Proteomic responses to hypoxia at different temperatures in the great scallop (Pecten maximus). PeerJ 2015, e871 (2015).
Mohamed, M. M. & Sloane, B. F. multifunctional enzymes in cancer. Nat. Rev. Cancer 2006 610 6, 764–775 (2006).
Pollutri, D. & Penzo, M. Ribosomal Protein L10: From Function to Dysfunction. Cells 2020, Vol. 9, Page 2503 9, 2503 (2020).
Guimaraes, J. C. & Zavolan, M. Patterns of ribosomal protein expression specify normal and malignant human cells. Genome Biol. 17, 1–13 (2016).
Mongin, A. A. & Orlov, S. N. MECHANISMS OF CELL VOLUME REGULATION. Physiol. Maint. 1, 130 (2009).
Sardini, A. et al. Cell volume regulation and swelling-activated chloride channels. Biochim. Biophys. Acta - Biomembr. 1618, 153–162 (2003).
De Zwaan, A., & Mathieu, M. Cellular Energy Metabolism in the Mytilidae: An Overview. in The Mussel Mytilus: Ecology, Physiology, Genetics and Culture (ed. Gosling, E.) 223–307 (Elsevier B.V., 1992).
DeBerardinis, R. J. & Chandel, N. S. We need to talk about the Warburg effect. Nat. Metab. 2020 22 2, 127–129 (2020).
Tseng, Y. C. et al. Regulation of lactate dehydrogenase in tilapia (Oreochromis mossambicus) gills during acclimation to salinity challenge. Zool. Stud. 47, 473–480 (2008).
Hochachka, P. W. Metabolic arrest. Intensive Care Med. 1986 123 12, 127–133 (1986).
Storey, K. B. & Storey, J. M. Metabolic rate depression and biochemical adaptation in anaerobiosis, hibernation and estivation. Q. Rev. Biol. 65, 145–174 (1990).
Hawkins, A. J. S. Relationships between the synthesis and breakdown of protein, dietary absorption and turnovers of nitrogen and carbon in the blue mussel, Mytilus edulis L. Oecologia 1985 661 66, 42–49 (1985).
Ramsden, D. B. et al. Human neuronal uncoupling proteins 4 and 5 (UCP4 and UCP5): structural properties, regulation, and physiological role in protection against oxidative stress and mitochondrial dysfunction. Brain Behav. 2, 468–478 (2012).
Papa, S. & Skulachev, V. P. Reactive oxygen species, mitochondria, apoptosis and aging. Detect. Mitochondrial Dis. 305–319 (1997) doi:10.1007/978-1-4615-6111-8_47.
Azzu, V. & Brand, M. D. The on-off switches of the mitochondrial uncoupling proteins. Trends Biochem. Sci. 35, 298–307 (2010).
Speed, N. & Blair, I. A. Cyclooxygenase- and lipoxygenase-mediated DNA damage. Cancer Metastasis Rev. 2011 303 30, 437–447 (2011).
Hove-Jensen, B. et al. Phosphoribosyl Diphosphate (PRPP): Biosynthesis, Enzymology, Utilization, and Metabolic Significance. Microbiol. Mol. Biol. Rev. 81, (2017).
Frances, A. & Cordelier, P. The Emerging Role of Cytidine Deaminase in Human Diseases: A New Opportunity for Therapy? Mol. Ther. 28, 357–366 (2020).
Zhao, X., Yu, H., Kong, L. & Li, Q. Transcriptomic Responses to Salinity Stress in the Pacific Oyster Crassostrea gigas. PLoS One 7, e46244 (2012).
Lockwood, B. L. & Somero, G. N. Transcriptomic responses to salinity stress in invasive and native blue mussels (genus Mytilus). Mol. Ecol. 20, 517–529 (2011).
Lehninger, A. L., Nelson, D. L., & Cox, M. M. Lehninger principles of biochemistry. (W.H. Freeman, 2005).
Koppenol, W. H., Bounds, P. L. & Dang, C. V. Otto Warburg’s contributions to current concepts of cancer metabolism. Nat. Rev. Cancer 2011 115 11, 325–337 (2011).
Zu, X. L. & Guppy, M. Cancer metabolism: facts, fantasy, and fiction. Biochem. Biophys. Res. Commun. 313, 459–465 (2004).
Su, M. A. et al. An Invertebrate Warburg Effect: A Shrimp Virus Achieves Successful Replication by Altering the Host Metabolome via the PI3K-Akt-mTOR Pathway. PLOS Pathog. 10, e1004196 (2014).
Lunt, S. Y. & Vander Heiden, M. G. Aerobic glycolysis: Meeting the metabolic requirements of cell proliferation. Annu. Rev. Cell Dev. Biol. 27, 441–464 (2011).
Leithner, K. PEPCK in cancer cell starvation. Oncoscience 2, 805 (2015).
Leithner, K., Hrzenjak, A. & Olschewski, H. Gluconeogenesis in cancer: Door wide open. Proc. Natl. Acad. Sci. U. S. A. 111, (2014).
Zhang, P. et al. Tumor suppressor p53 cooperates with SIRT6 to regulate gluconeogenesis by promoting FoxO1 nuclear exclusion. Proc. Natl. Acad. Sci. U. S. A. 111, 10684–10689 (2014).
Evan, G. I. & Vousden, K. H. Proliferation, cell cycle and apoptosis in cancer. Nat. 2001 4116835 411, 342–348 (2001).

Tables 2 to 5 are available in the Supplementary Files section.

No competing interests reported.

Download PDF

Version 1

posted

You are reading this latest preprint version

Transcriptional Responses to Long-term Salinity Stress and Acclimation in the Black Mussel, Mytilus Galloprovincialis

Status:

Version 1

Abstract

Figures

Introduction

Materials And Methods

Results

Discussion

Declarations

References

Tables

Additional Declarations

Supplementary Files

Status:

Version 1