Hypoosmotic stress induced tissue-specific immune responses of yellowfin seabream (Acanthopagrus latus) revealed by transcriptomic analysis

Journal Pre-proof Hypoosmotic stress induced tissue-specific immune responses of yellowfin seabream (Acanthopagrus latus) revealed by transcriptomic analysis Genmei Lin, Min Zheng, Dong Gao, Shizhu Li, Wenyu Fang, Jing Huang, Jingui Xie, Jingxiong Liu, Yijing Liu, Zhaohong Li, Jianguo Lu PII:

S1050-4648(20)30114-5

DOI:

https://doi.org/10.1016/j.fsi.2020.02.028

Reference:

YFSIM 6827

To appear in:

Fish and Shellfish Immunology

Received Date: 9 October 2019 Revised Date:

12 February 2020

Accepted Date: 14 February 2020

Please cite this article as: Lin G, Zheng M, Gao D, Li S, Fang W, Huang J, Xie J, Liu J, Liu Y, Li Z, Lu J, Hypoosmotic stress induced tissue-specific immune responses of yellowfin seabream (Acanthopagrus latus) revealed by transcriptomic analysis, Fish and Shellfish Immunology (2020), doi: https://doi.org/10.1016/j.fsi.2020.02.028. This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2020 Published by Elsevier Ltd.

CRediT author statement Jianguo Lu: Conceptualization, Methodology, Resources, Writing - Review & Editing, Supervision, Project Administration, Funding Acquisition. Genmei Lin: Methodology, Validation, Formal Analysis, Investigation, Data Curation, Writing Original Draft, Writing - Review & Editing, Visualization. Min Zheng: Formal Analysis, Investigation, Data Curation, Writing - Review & Editing. Shizhu Li: Methodology. Dong Gao: Methodology. Wenyu Fang: Validation. Jing Huang: Investigation. Jingui Xie: Investigation. Jingxiong Liu: Investigation. Yijing Liu: Investigation. Zhaohong Li: Investigation.

Graphical Abstract:

1

Hypoosmotic stress induced tissue-specific immune responses of

2

yellowfin seabream (Acanthopagrus latus) revealed by

3

transcriptomic analysis

4 5

Genmei Lina,1, Min Zhenga,1, Dong Gaoa, Shizhu Lia, Wenyu Fanga, Jing Huanga,

6

Jingui Xiea, Jingxiong Liua, Yijing Liua, Zhaohong Lia, Jianguo Lua,b,*

7

a

School of Marine Sciences, Sun Yat-sen University, Zhuhai, 519082, China

8

b

Southern Marine Science and Engineering Guangdong Laboratory (Zhuhai), Zhuhai

9

519080, China

10

1

Authors contributed equally to this work.

11

*

Corresponding author: [email protected]

12

1

13 14

Abstract Salinity is a limiting factor for many marine organisms, including fishes. The

15

shift in the ambient salinity can cause osmotic stress and arouse immune responses in

16

fish. In this study, yellowfin seabream (Acanthopagrus latus), a euryhaline marine

17

teleost, was used to investigate immune responses of different tissues (gill, liver, and

18

muscle) under hypoosmotic stress. Comparative transcriptomic and physiological

19

analyses of three tissues were conducted after fish exposed to the fresh water (FW,

20

salinity = 0 ppt), low-saline water (LW, salinity = 3 ppt), and brackish water (BW,

21

salinity = 6 ppt) for 8 days. The results showed that hypoosmotic stress dramatically

22

altered the gene expression of three tissues in yellowfin seabream; The investigation

23

of differentially expressed genes (DEGs) related to osmoregulation and immune

24

response indicated that T cell-mediate immunity pathways were essential to tackle

25

such stress. In terms of tissues, gill was found to be the most sensitive tissue under

26

hypoosmotic stress by enhancing of Na+K+-ATPase activity and preventing the loss of

27

Na+ and K+; Liver, on the other hand, was under the most sever oxidative stress

28

indicated by the fluctuation of SOD, CAT activities and the MDA content; In contrast,

29

muscle had the least osmoregulation and immune related response. We also identified

30

several potential candidate genes, which may serve as gene indicators to identify the

31

stressor. Overall, this study provides preliminary mechanistic insights into

32

hypoosmotic stress adaption of aquatic organism.

33

Keywords: Yellowfin seabream; Hypoosmotic stress; Osmoregulation; Immune

34

responses; RNA-seq analysis 2

35 36

1. Introduction Salinity is a major abiotic factor in the environment for all aquatic organisms (Su

37

et al., 2019). A change in the salinity of habitat water can lead to osmotic stress and

38

influence hormonal levels, metabolic and energetic parameters and osmoregulation or

39

immune response, etc. in fish (Boeuf and Payan, 2001; Bowden, 2008;

40

Vargas-Chacoff et al., 2015). One of the most important regulatory systems to tackle

41

such stress is immune system. Stress induced by salinity changes can agitate the

42

generation of reactive oxygen species (ROS) that may compromise the efficacy of

43

immune system and cause oxidative stress (Tort, 2011). Marine fish can also

44

recognize and respond to salinity variation by altering the expression of genes

45

specifically required for osmoregulation (Su et al., 2019).

46

To identify such abiotic stress, researchers reported a series of innovative

47

analytical approaches, including flow cytometry, automatic blood cell analysis,

48

micro-Raman spectroscopy and 1H NMR and immunohistochemical assays to assess

49

the health status of aquaculture species under the stress condition such as diseases,

50

environmental changes and presence of pollutants (Acar et al., 2018; Cappello et al.,

51

2015; Fazio et al., 2019; Parrino et al., 2019a, 2019b). For example, analyzing the

52

whole-blood samples of fish by flow cytometry could tell that the increasing water

53

salinity induced erythropoiesis in teleost fish species (Parrino et al., 2018). If the

54

underpinning genes and related mechanisms are concerned, transcriptomics is the

55

most researched omics technology for such purpose.

56

Most of studies related to osmotic stress in fish have been focus on the 3

57

hyperosmotic stress (Birrer et al., 2012; Cuesta et al., 2005; Jiang et al., 2008; Marc et

58

al., 1995; Taylor et al., 2007), in which they stated that acute or chronic salinity

59

changes could inhibit the expression of some immune-related genes and arouse both

60

cellular and humoral immune activities. Nonetheless, hypoosmotic stress could also

61

modulate the immune homeostasis of fish (Schmitz et al., 2016).

62

As euryhalines, the fish must be able to cope with salt depletion and water gain

63

when it is in fresh water. It was reported that when gilthead seabream (Sparus aurata)

64

was acclimated to water with the low salinity of 6 ppt, it exhibited lower peroxidase

65

and alternative complement activity in plasma compared to those acclimated to water

66

with higher salinity (12 ppt and 38 ppt) (Cuesta et al., 2005). Another study stated that

67

broad-nosed pipefish (Syngnathus typhle) at low salinity were unable to mount

68

specific immune response components upon bacterial infection (Birrer et al., 2012).

69

However, studies evaluating fish acclimation to hypoosmotic environment and the

70

underline mechanism of such stressors on the immune capacity remain largely

71

unexplored.

72

Yellowfin seabream, Acanthopagrus latus (Houttuyn), is an economically and

73

ecologically important marine fish which is broadly distributed throughout the

74

Indo-West Pacific area including southern Japan, southeastern China, Taiwan, Persian

75

Gulf, and Australia (Abou-Seedo et al., 2003; Nakabo, 2002; Xia et al., 2006, 2008).

76

Owing to its excellent quality and good taste, it has been successfully farmed in

77

Southeast Asia since the 1980s (Hu et al., 2007). The salinity level for offshore

78

culture in south China has been reported within the range of 4 to10 ppt for this species 4

79

(Hu et al., 2007; Wu et al., 2007a, 2007b). Inevitably, yellowfin seabream could

80

encounter hypoosmotic stress when moving from coastal waters to offshore brackish

81

water culture environment. Earlier studies related to salinity adaption of yellowfin

82

seabream had been focused on their morphological (Movahedinia et al., 2009) and

83

physiological (Farshadian et al., 2018) responses. However, the underline adaption

84

molecular mechanism especially for hypoosmotic stress is unclear.

85

In this study, we employed the next-generation sequencing (NGS) technology to

86

sequence and de novo assembled the gill, liver, and muscle transcriptome of yellowfin

87

seabream in the brackish water (BW, salinity = 6 ppt), the low-saline water (LW,

88

salinity = 3 ppt), and the fresh water (FW, salinity = 0 ppt) environments. The aim of

89

this study is to identify immune system related transcriptomic differences associated

90

with FW- and LW-acclimated tissues within a single species along with the related

91

metabolic pathways. To our best knowledge, this is the first NGS study that had

92

provided a differential transcriptome-wide perspective among tissues of a euryhaline

93

fish acclimated to FW and LW. Also, these candidate immune related genes may be

94

used as molecular biomarkers to diagnose stress in fish, improve stress tolerance and

95

reduce diseases of fish in the activities of aquaculture and fisheries stock

96

enhancement.

97

2. Material and methods

98

2.1 Fish maintenance and salinity challenge experiment

99 100

All protocols were carried out in accordance with the Ethic Committee on the Care and Use of Laboratory Animals of Sun Yat-sen University, Guangzhou, China. 5

101

Yellowfin seabream with average body length of 11.2 ± 0.9 cm and body weight

102

of 36.93 ± 10.57 g were obtained from a local fish farm (Zhuhai, China) at a salinity

103

of 6.06 ± 0.47 ppt, Then the fish were transferred to our lab and kept in the tank filled

104

with brackish water (BW) in the salinity of 5.98 ± 0.09 ppt for two weeks for adaption

105

to the new condition before exposing them to other salinity variations. Fish were

106

maintained at 26.5 ± 0.5°C under a 14h/10h light-dark photoperiod with automated

107

water circulation and gentle aeration. During this period, fish were fed with red worm

108

twice daily. No mortality was recorded during the entire acclimation period.

109

Three salinity including ~ 6, 3 and 0 ppt were used in this study in triplicate. Fish

110

were fasted 24h before salinity challenge. After acclimation, fish were suddenly

111

exposed to fresh water (FW, the salinity was adjusted to 0.10 ± 0.03 ppt) or low-saline

112

water (LW, the salinity was adjusted to 2.94 ± 0.40 ppt) respectively (n = 12 per tank).

113

One group remained in the brackish water (n = 12) as a control. Each tank was filled

114

with 54L water under the culture conditions mentioned above. After 8 days (Xia et al.,

115

2013), the survival rates of yellowfin seabream were 100%, 100% and 58% in BW,

116

LW, and FW group respectively.

117

Three fishes per tank were anesthetized with tricaine methane sulfonate (MS 222,

118

200 mg/L) and sampled immediately. The tissues (gill, liver, and muscle) were

119

collected and stored at -80°C till use.

120

2.2 Physiological parameter analyses

121 122

Gill, liver, and muscle tissues from the same individual from each group were collected for physiological analyses. Tissues were homogenized in 1.0 mL deionized 6

123

water for Na+ concentration, K+ concentration, and Na+K+-ATPase (NKA) activity

124

measurements or normal saline for superoxide dismutase (SOD) activity, catalase

125

(CAT) activity, and the malonic dialdehyde (MDA) contents measurements on ice.

126

After centrifugation, the supernatant was transferred to a new tube and used for

127

physiological analyses. Commercially kits were used to test Na+ concentration, K+

128

concentration, NKA activity, SOD activity, CAT activity and the MDA contents of

129

three tissues following the manufacturer's protocol (Nanjing Jiancheng

130

Bioengineering Institute, China). The total protein content was analyzed following the

131

iodine-starch colorimetric method.

132

2.3 RNA extraction, library construction, and transcriptome sequencing

133

Total RNA was extracted using TRIzol Reagent (Invitrogen, USA) following

134

manufacturer’s instructions. Briefly, the tissues were homogenized using bead-beating

135

homogenizer (TOMOS, China) and resuspended in 1 mL TRIzol solution. The

136

mixture was shaken, incubated at room temperature for 5 min and centrifuged at 4 °C.

137

Then the supernatant was removed to a new tube and mixed with 200 µL chloroform.

138

The mixture was vortexed, incubated at room temperature for 5 min, and centrifuged

139

for 10 min at what 10,000 × g at 4 °C. The supernatant was mixed with 500 µL

140

isopropanol and stored at -20 ℃ for 1 h and then kept at room temperature for 10 min.

141

The RNA was precipitated and washed twice with 1 mL 70% RNase-free ethanol,

142

then air dried, and resuspended in 100 µL DEPC (diethyl pyrocarbonate)-treated

143

water.

144

The concentration and purity of extracted RNA was checked using a NanoDrop 7

145

2000 Spectrophotometer (Thermo Fisher Scientific, USA). The quality of RNA was

146

checked through electrophoresis in 1% agarose gel to confirm there was no

147

degradation and genomic DNA contamination. The integrity of RNA was assessed

148

using the RNA Nano 6000 Assay Kit of the Agilent Bioanalyzer 2100 system (Agilent

149

Technologies, USA).

150

RNA with RIN value above 8.0 was used for library construction. The library

151

was generated using TruseqTM RNA Sample Prep Kit (Illumina, USA) following

152

manufacturer’s instructions. Briefly, mRNA was purified from the total using

153

oligo(dT)-attached magnetic beads, which was then fragmented using fragmentation

154

buffer. The first-strand cDNA was synthesized using random hexamer primers and

155

reverse transcriptase. The second-strand cDNA was synthesized subsequently using

156

RNase H and DNA polymerase I. The double stranded cDNA was polished blunt

157

ended with end repair mix. After dA addition at the 3’ ends of cDNA fragments,

158

adaptors with hairpin loop structure were added. The cDNA was enriched and

159

amplified through PCR. The PCR product was purified, and its quality was assessed

160

on Bioanalyzer 2100 system (Agilent, USA). The libraries were pair ends sequenced

161

on Illumina HiseqTM platform (Shanghai Majorbio Bio-pharm Technology Co., Ltd,

162

China).

163

2.4 Bioinformatic analyses

164

Raw reads were first trimmed and quality-filtered using Trimmomatic v0.38

165

(Bolger et al., 2014). GC content, Q20, Q30, and sequence duplication level were

166

calculated through FastQC v0.11.8. Transcriptome was assembled using Trinity v2.8.4 8

167

(Grabherr et al., 2011) with default parameters. The Trinity assembly was optimized

168

using CD-HIT v4.6.8 (Li and Godzik, 2006) and TransRate v1.0.3 (Smith-Unna et al.,

169

2016) and then evaluated using BUSCO v3 (Simão et al., 2015). Databases including

170

NR (NCBI non-redundant protein sequences), GO (Gene ontology), KEGG (Kyoto

171

encyclopedia of genes and genomes ortholog database), KOG (Clusters of eukaryotic

172

orthologous groups of proteins), Pfam (Protein family), and Swiss-Prot (A manually

173

annotated and reviewed protein sequence database) were applied to gene annotation.

174

Transcript abundance was quantified and normalized into Transcripts Per Million

175

reads (TPM) with RSEM v1.2.15 (Li and Dewey, 2011). Differential expression of

176

genes between libraries was analyzed using the DESeq2 R package (Love et al., 2014).

177

The threshold for differentially expressed genes (DEGs) was p-adjust value < 0.005

178

and log2|foldchange| > 1. For functional annotation and classification, the DEGs were

179

enriched into GO terms and KEGG pathways using Goatools (Klopfenstein et al.,

180

2018) and KOBAS (Mao et al., 2005), with p-adjust value < 0.05 considered

181

statistically significant.

182

The six physiological parameters (including Na+, K+ concentration, NKA activity

183

SOD, CAT activity, and the MDA contents) and total transcript abundance of each

184

KEGG pathways under “immune system” category were used for correlation analysis

185

by R using rcor.test (Rizopoulos, 2006) and with “qvalue” package (Storey et al.,

186

2019) for false positives test.

187

2.5 qRT-PCR validation

188

To assess the expression level of randomly selected DEGs, extracted RNA used 9

189

for transcriptome sequencing was also used to perform qRT-PCR. Briefly, total RNA

190

was reverse transcribed into cDNA using PrimeScript™ RT reagent Kit with gDNA

191

Eraser (TaKaRa, China). qRT-PCR was carried out on a LightCycler®️ 480 II System

192

(Roche Diagnostics, Switzerland) using TBGreen®️ Premix Ex TaqTM II (TaKaRa,

193

China). The primer sequences for the target genes are listed in SI Table S1. The cycle

194

threshold (Ct) values of target genes were normalized to housekeeping gene β-actin

195

expression level as internal control using the 2-△△Ct method. All qRT-PCR experiments

196

were performed in triplicate.

197

2.6 Statistical analyses

198

The significant differences of physiological parameters among each group were

199

performed by one-way analysis of variance (ANOVA) with SPSS Version 21.0.

200

Values with p < 0.05 were considered significant. All values are expressed as the

201

mean ± standard deviation.

202

3. Results

203

3.1 Physiological analyses of the gills, livers and muscles

204

3.1.1 Na+ concentration, K+ concentration, and Na+K+-ATPase activity

205

The average of Na+ and K+ concentration in gills, livers, and muscles were

206

measured in this study, and the values ranged from 2.40 ± 0.55 to 4.20 ± 0.50

207

mmol/gprot for concentration of Na+, and 0.25 ± 0.07 to 1.08 ± 0.32 mmol/gprot for

208

concentration of K+ in three tissues. Gills harbored the highest Na+ and K+

209

concentration, whereas, muscles and livers had the lowest Na+ and K+ concentration,

210

respectively. However, no significant differences of Na+ and K+ concentration were 10

211

detected among three treatments in three tissues (Figure 1A; Figure 1B; SI Table

212

S2).

213

Na+K+-ATPase (NKA), a kind of protease existing on the membrane of cells and

214

organelles, plays a pivotal role in substances transportation, energy transformation,

215

and information transmission. The activities of NKA in three tissues were further

216

investigated. In gills, the values of NKA activity in FW group (3.52 ± 0.83 U/mgprot,

217

p-value < 0.01) were found significantly increased compared with that in LW group

218

(2.04 ± 0.53 U/mgprot) and BW group (2.31 ± 0.60 U/mgprot) Nonetheless, no

219

significant differences were found for the average value of NKA activity among three

220

treatments in livers and muscles (liver: 0.19 ± 0.03 U/mgprot; muscle: 1.00 ± 0.20

221

U/mgprot) (Figure 1C).

222

3.1.2 SOD activity, CAT activity and the MDA contents

223

As shown by the results (Figure 1D; Figure 1E), the antioxidant capability of

224

yellowfin seabream was significantly compromised upon FW exposure in gills, as

225

indicated by the significant inhibition of SOD activity (71.2% of SOD activity of BW

226

exposure group with p-value < 0.05) and reduction of CAT activity (35.8% of CAT

227

activity of BW exposure group) in the FW group. While, in LW group, SOD and CAT

228

activities were enhanced 1.46-fold and 1.30-fold compared with those in BW group.

229

Furthermore, the MDA levels significantly increased in FW (1.75-fold) and LW

230

(1.30-fold) groups compared with those in BW group (Figure 1F).

231

The SOD, CAT activities, and the MDA contents in livers had more obvious

232

responses (p < 0.01) than those in gills across treatments, yet with similar regulation 11

233

patterns (Figure 1D; Figure 1E; Figure 1F). It demonstrated the importance of liver

234

in osmoregulation of yellowfin seabream. Whereas, in muscles, only SOD activity

235

was significantly depressed about 40.5% and 4.4% in FW and LW group compared

236

with those in BW group. In contrast, CAT activity and the MDA contents showed no

237

significant changes across three groups (Figure 1D; Figure 1E; Figure 1F).

A

B

C

6

2

6

Liver

Mu scle

E

NKA (U/mgprot)

F

40

* ** **

10

20

0

0 Liver

MDA (nmol/mgprot)

CAT (U/mgprot)

** **

20

60

Mu scle

brackish water (6 ppt)

**

* **

30

Gill

Liver

Mu scle

50

*

SOD (U/mgprot)

Mu scle

*

240 241 242 243 244

Liver

100

Gill

239

Gill

B

80

2

0

0 Gill

50 40

0.5

4

40 30 20

*

0

D

1

**

2

1.5

**

4

*

K+ (mmol/gpot)

B

** **

Na+ (mmol/gprot)

238

10 0

Gill

Liver

low-saline water (3 ppt)

Mu scle

Gill

Liver

Muscle

fresh water (0 ppt)

Figure 1 Physiological parameters of gills, livers and muscles in yellowfin seabream exposed to fresh, low-saline and brackish water for 8 days. (A) Na+ concentrations; (B) K+ concentrations; (C) Na+K+-ATPase activities; (D) SOD activities; (E) CAT activities; and (F) The MDA contents. * p < 0.05 and ** p < 0.01 indicate significant difference between each group.

245

3.2 Transcriptomic analyses of the gills, livers and muscles

246

3.2.1 Transcriptome de novo assembly and functional annotation

247

The results revealed that the salinity change induced dramatic alterations in

248

global gene transcription profiles in three tissues of yellowfin seabream. It produced

249

approximately a total of 77.69 Gbp, 70.73 Gbp, and 86.32 Gbp clean sequence data

250

for the gill, liver, and muscle tissue samples, respectively, with an average GC

251

percentage of 50.03%. The average of Q20 and Q30 were 97.52% and 92.74%, 12

252

respectively (SI Table S3). Clean reads were assembled into 102,215 unigenes

253

(hereafter genes). The lengths ranged from 201 bp to 26,594 bp, with an N50 length

254

of 2,089 bp and an average length of 1,270 bp (SI Figure S1).

255

Of these genes, 35,697 (34.92%) were annotated by referring to the currently

256

available databases including NR, GO, KO, Swiss-Prot, Pfam, and KOG (Table 1);

257

We utilized GO assignment programs for functional categorization of those annotated

258

genes, and they were able to map 12,329 genes (12.06%) to GO terms. In many cases,

259

multiple terms were assigned to the same transcript. This resulted in assignments of

260

7,347 genes to cellular components (CC), 6,364 to biological processes (BP), and

261

10,227 to molecular functions (MF) (SI Figure S2). A total of 22,987 genes (22.49%)

262

involved in 43 KEGG secondary metabolic pathways were detected in our

263

transcriptome database (SI Figure S3). In addition, 12,701 genes (12.43%) were

264

divided into 24 KOG clusters (SI Figure S4).

265 266

267 268 269 270 271

Table 1 Summary of annotated unigenes of yellowfin seabream. Database

Number of annotated unigenes

Percentage of annotated unigenes (%)

NR

33,673

32.94

Swiss-Prot

28,554

27.94

Pfam

25,335

24.79

KOG

12,701

12.43

GO

12,329

12.06

KEGG

22,987

22.49

Annotated in at least one database

35,697

34.92

Total

102,215

100.00

Note: NR, NCBI non-redundant protein sequences; Swiss-Prot, a manually annotated and reviewed protein sequence database; Pfam, Protein family; KOG, Clusters of eukaryotic orthologous groups of proteins; GO, Gene ontology; KEGG, Kyoto encyclopedia of genes and genomes ortholog database.

In order to identify transcripts conserved in the 3 tissues (gill, liver, and muscle)

13

272

and those unique to each tissue under hypoosmotic stress, a comparative analysis

273

among the yellowfin seabream transcriptomes was conducted (Figure 2).

274

In gill, a total of 1,347 genes (up-regulated DEGs: 518, 38.46%, down-regulated

275

DEGs: 829, 61.54%) were identified as DEGs in FW group compared to BW group

276

(Figure 2A), and significantly enriched into 53 GO terms (SI Table S4) and 33

277

KEGG pathways (SI Table S5). In contrast, only 10 and 17 DEGs were found up and

278

down regulated in LW group compared to BW group respectively (Figure 2D). A total

279

of 27 DEGs were enriched in 7 KEGG pathways (SI Table S5). However, no GO

280

terms were enriched in this case (SI Table S4).

281

In liver, compared with BW, 1,116 genes were differentially expressed in FW.

282

Among them, 782 genes (70.07%) were up-regulated, whereas, 334 genes (29.93%)

283

were down-regulated (Figure 2B). These DEGs were further enriched into 13 KEGG

284

pathways, while no significantly enriched GO terms were found (SI Table S6; SI

285

Table S7). For comparison, a total of 525 DEGs (up-regulated DEGs: 171, 32.57%

286

and down-regulated DEGs: 354, 67.43%) were found in LW group compared with

287

those in BW group (Figure 2E), which were identified enriched into 19 GO terms and

288

4 KEGG pathways (SI Table S6; SI Table S7).

289

In muscle, a total of 1,584 DEGs (728, 45.96% up-regulated and 856, 54.04%

290

down-regulated) were found in FW group compared with those in BW group (Figure

291

1C), which were further significantly enriched in 7 GO terms and 8 KEGG pathways

292

(SI Table S8; SI Table S9). In LW group, on the other hand, fewer DEGs were found

293

compared with BW. Namely, only 33 (19, 57.58% up-regulated and 14, 42.42% 14

294

down-regulated) DEGs were identified (Figure 1F) and enriched in 2 KEGG

295

pathways but not in any GO terms (SI Table S8; SI Table S9).

296

297 298 299 300 301 302 303 304 305 306 307 308

Figure 2 Volcano plots for gill, liver and muscle gene libraries of fresh water and brackish water, or low-saline water and brackish water of yellowfin seabream by showing variance in gene expression with respect to fold change (FC) and significance (p-value). (A) Volcano plots for gill gene libraries of fresh water and brackish water; (B) Volcano plots for liver gene libraries of fresh water and brackish water; (C) Volcano plots for muscle gene libraries of fresh water and brackish water; (D) Volcano plots for gill gene libraries of low-saline water and brackish water; (E) Volcano plots for liver gene libraries of low-saline water and brackish water; and (F) Volcano plots for muscle gene libraries of low-saline water and brackish water. Each dot represents an individual gene: Red dots represent the down-regulated DEGs, blue dots represent the up-regulated DEGs, and grey dots represent not differentially expressed genes.

3.2.2 Analysis of DEGs for osmoregulation in three tissues of yellowfin seabream

309

Specifically, DEGs related to osmoregulation in three tissues of yellowfin

310

seabream were further identified in three tissues and mainly classified into three GO

311

terms, including regulation of ion transmembrane transport (GO:0034765) (Figure

312

3A), integral component of membrane (GO:0016021) (Figure 3B), and ATP binding

313

(GO:0005524) (Figure 3C). Description of selected DEGs was listed in SI Table S11.

314

In teleost fish, gill is one of the major tissues to maintain osmotic balance by

315

coordinated water and ion transport (Evans, 2008; Laurent and Dunel, 1980). As 15

316

shown in our data, gill harbored more osmoregulation related DEGs than liver and

317

muscle (DEGs: 23 in gill, 21 in liver, and 17 in muscle) (Figure 3). DEGs including

318

kcnj15 (ATP-sensitive inward rectifier potassium channel 15) and atp (Na+K+-ATPase

319

alpha subunit) showed significantly up-regulation, while clic4 (chloride intracellular

320

channel protein 4) and ncp4 (sodium channel protein type 4 subunit alpha-like) were

321

significantly down-regulated in FW group compared with BW group (Figure 3).

322

Strikingly, significant up-regulation of hsp70 (heat shock protein 70) multigene family

323

were identified in gill (about 17.89 and 22.25 folds) in FW group compared with that

324

in BW group (Figure 3).

325

Although liver and muscle are not the main osmo-regulatory tissues, several

326

classical ion transporters were also found in our DEG list. To our surprise, the

327

expression of slc22a16 (solute carrier family 22 member 16) was found about 92-fold

328

up-regulated in the liver of FW group compared with that in BW group (Figure 3B).

329

Another DEG worthwhile to mention is cret1l (sodium- and chloride-dependent

330

creatine transporter 1-like) which performed similarly to slc22a16 and showed

331

approximately 61-fold over-expression under hypoosmotic stress (FW group) (Figure

332

3B). In addition, hsp70 genes were found significantly over-expressed in liver (35.67

333

and 47.21-fold) in FW group compared with those in BW group as well (Figure 3C).

334

In muscle, multiple genes were found significantly down-regulated by approximately

335

6- to 9-fold, such as ncp4, kcna7 (potassium voltage-gated channel subfamily A

336

member 7-like), aqp7 (aquaporin-7), and grp78 (78 kDa glucose-regulated protein) in

337

FW compared with those in BW group. Whereas, slc22a16, cret1l, ntat (large neutral 16

338

amino acids transporter small subunit 2-like isoform X2), and nslc3 (sodium-coupled

339

neutral amino acid transporter 3), were found significantly up-regulated by about 6- to

340

11-fold in FW group compared with their counterparts in BW group.

341

It is noteworthy that the DEGs associated with osmoregulation in LW group

342

showed the same regulation pattern but with lower fold change values compared with

343

those identified in FW group.

344

345 346 347 348 349 350 351 352 353

Figure 3 Heatmap of differentially expressed genes (DEGs) annotated into osmoregulation-related GO terms in gill, liver, and muscle of yellowfin seabream exposed to fresh, low-saline and brackish water. (A) DEGs annotated into regulation of ion transmembrane transport (BP, GO:0034765); (B) DEGs annotated into integral component of membrane (CC, GO:0016021); and (C) DEGs annotated into ATP binding (MF, GO:0005524). BP: biological process; CC: cellular component; MF: molecular function. Note: Data shown had the absolute fold change > 2 at p <0.05. Transcript abundance was quantified and normalized into Transcripts Per Million reads (TPM) and (TPM+1) was calculated as logs base ten.

354

3.2.3 Analysis of DEGs for immune responses in three tissues of yellowfin

355

seabream

356

Furthermore, immune-related DEGs were identified based on “immune system” 17

357

category in KEGG pathway database. Overall, immune related DEGs had the same

358

regulation pattern in FW and LW groups compared with those in BW group, and more

359

profound fold changes were found in FW group (fold change: 4.17 to 29.13) than LW

360

group (fold change: 1.00 to 2.04) in three tissues compared with BW group, yet with

361

only one exception, that is to say, c8g (complement component C8 gamma chain)

362

strikingly overexpressed in LW (10.45-fold) rather than in FW group (1.00-fold)

363

(Figure 4).

364

A series of immune-related DEGs play crucial roles in response of hypoosmotic

365

stress in yellowfin seabream. For example, interleukin encoding DEGs, il12p40

366

(interleukin 12p40), il17af (interleukin 17 A/F), and il17ra1 (interleukin 17 receptor

367

A1) were identified down-regulated (9.09-fold, 11.63-fold, and 4.35-fold inhibition)

368

in gills when exposed to FW group (Figure 4), while il21 (interleukin-21 receptor-like)

369

and il20 (interleukin-20 receptor subunit alpha-like) were found over-expressed to

370

8.24-fold and 4.93-fold in livers under FW exposure compared with those in BW

371

group (Figure 4). Multiple isoforms of complement components which were

372

significantly enriched to complement and coagulation cascades pathway were

373

identified up-regulated in this study as well, such as c8g in gill and muscle, f8l

374

(coagulation factor VIII-like) and f9l (coagulation factor IX-like) in liver (Figure 4).

375

The immune-related KEGG pathway enrichment analysis was performed based

376

on the loading stress-responsive genes to identify the significantly stress-related

377

pathways. Overall, the majority of immune-related KEGG pathways were enriched in

378

gill, followed by liver and then muscle. In gill, KEGG pathways including Th1 and 18

379

Th2 cell differentiation (ko04658), IL-17 signaling pathway (ko04657), T cell

380

receptor signaling pathway (ko04660), and natural killer cell mediated cytotoxicity

381

(ko04650) were significantly enriched. Jak-STAT signaling pathway (ko04630) and

382

NOD-like receptor signaling pathway (ko04621) were significantly enriched in liver,

383

while B cell receptor signaling pathway (ko04662) was enriched in muscle.

384

Furthermore, Th17 cell differentiation (ko04659) pathway was significantly enriched

385

in both gill and liver, and antigen processing and presentation (ko04612) pathways

386

were significantly enriched in both gill and muscle (Figure 4). Besides, several

387

immune related DEGs, including btk (tyrosine-protein kinase BTK isoform X1) and

388

blnk (B-cell linker protein) in gill, cd4 (T cell surface glycoprotein CD4) and ccl4l

389

(C-C motif chemokine 4-like) in liver, as well as gabarapl1 (gamma-aminobutyric

390

acid receptor-associated protein-like), samd9l (sterile alpha motif domain-containing

391

protein 9-like), and c8g in muscle were found under or over expressed with higher

392

fold change value (absolute value > 4) under salinity stress, however they didn’t

393

evoke any KEGG pathway regulation significantly (Figure 4).

394

19

fresh water vs. brackish water

Gill

low-saline water vs. brackish water

-10

0

10

20

-10

0

10

20

-20

395 396 397 398 399 400 401 402 403 404

-10

0

10

TRINITY_DN52172_c0_ g1 TRINITY_DN60540_c1_ g1 TRINITY_DN55961_c0_ g1 TRINITY_DN50631_c0_ g1 TRINITY_DN58309_c0_ g1 TRINITY_DN53214_c1_ g1 TRINITY_DN52161_c0_ g1 TRINITY_DN54733_c0_ g1 TRINITY_DN38306_c0_ g1 TRINITY_DN53985_c0_ g1 TRINITY_DN61500_c2_ g1 TRINITY_DN55790_c0_ g1 TRINITY_DN54724_c0_ g1

ahr1 zfcp gabarapl1 txnip il21 il20 stat5b cd4 ccl4l epo f8l c1r f9l

TRINITY_DN55961_c0_ g1 TRINITY_DN54706_c2_ g1 TRINITY_DN60879_c0_ g1 TRINITY_DN51818_c0_ g1 TRINITY_DN42389_c0_ g1 TRINITY_DN62310_c12_g1 TRINITY_DN57286_c0_ g1 TRINITY_DN54030_c0_ g1 TRINITY_DN51010_c1_ g1 TRINITY_DN44919_c0_ g1

gabarapl1 samd9l c8g nfatc nfatc mhca cd4 lgmn hsp70a4 hsp70a4l

A B C D E F G

30

Muscle

-20

il12p40 mcpt1l gzmbl gzm il17af ccl1l il17ra1 m2l crlf3 ccl17l m2l cd3e c8g btk blnk rt1-ab ht1-aqa rt1-ab hsp70 hsp70

30

Liver

-20

TRINITY_DN53015_c3_ g1 TRINITY_DN59006_c1_ g2 TRINITY_DN59006_c1_ g1 TRINITY_DN22548_c0_ g1 TRINITY_DN42824_c0_ g1 TRINITY_DN41743_c1_ g1 TRINITY_DN59828_c4_ g1 TRINITY_DN57391_c0_ g1 TRINITY_DN54733_c0_ g1 TRINITY_DN52530_c1_ g1 TRINITY_DN50208_c0_ g1 TRINITY_DN50899_c1_ g1 TRINITY_DN57870_c0_ g1 TRINITY_DN57311_c1_ g1 TRINITY_DN56925_c0_ g1 TRINITY_DN58789_c2_ g1 TRINITY_DN45722_c0_ g1 TRINITY_DN40780_c0_ g1 TRINITY_DN52028_c4_ g2 TRINITY_DN51512_c0_ g1

20

C F H I

30

Fold Change

A J

DEGs

Pathways

Figure 4 Fold changes of representative immune-related differentially expressed genes (DEGs) in gill, liver, and muscle and their significantly enriched KEGG pathways of yellowfin seabream exposed to fresh, low-saline and brackish water. KGGG pathways include, A: Antigen processing and presentation; B: Th1 and Th2 cell differentiation; C: Th17 cell differentiation; D: IL-17 signaling pathway; E: T cell receptor signaling pathway; F: Complement and coagulation cascades; G: Natural killer cell mediated cytotoxicity; H: Jak-STAT signaling pathway; I: NOD-like receptor signaling pathway; and J: B cell receptor signaling pathway. Note: Data shown had the absolute fold change > 4 at p <0.05.

3.2.4 Validation of RNA-seq by RT-qPCR

405

Several members of immune-related genes in three tissues were randomly

406

selected for further RT-qPCR validation, respectively. Fold changes represented the

407

ratio of relative gene expression levels between FW/LW and BW treatment (Figure 5;

408

SI Table S10). As an internal gene, the beta-actin was included. Under the given

409

experimental conditions, the exposure to three salinity concentration showed no effect

410

on the internal gene (data not shown). Generally, the gene expression levels and trend 20

411

of regulation (increase/decrease) detected with RNA-seq and qualification based on

412

RT-qPCR agreed with each other.

413 30

Gill

Liver

*

Fold Change

20

Muscle

10

0

rt1-ab htl-aqa

il17af

gzm

mcpt1l il12p40

ccl4l

hsp70a4l stat5b

il21

txnip gabarapl1 ahr1

nfatc c8g

cd4

samd9l gabarapl1

-10

-20 * *

-30

*

* *

fresh water vs. brackish water

low-saline water vs. brackish water

414 415 416 417 418 419 420 421 422 423

Figure 5 Validation of gill, liver and muscle tissue transcriptome results by RT-qPCR using 18 randomly selected differentially expressed genes upon exposure to yellowfin seabream exposed to fresh, low-saline and brackish water for 8 days. The RT-qPCR fold changes are relative to brackish water samples and are normalized by changes in beta-actin values. The mean of three relative quantities of biological replicates were used in a two-tailed Student's t-test with a 95% confidence level (p < 0.05) to determine the gene expression significance.

3.3 Correlation analysis between physiological parameters and immune system related KEGG pathways Correlation analysis demonstrated that immune system related KEGG pathways

424

regulation correlated with osmoregulation and antioxidant activities. However, the

425

correlation patterns varied across three tissues. In gill, six immune-related pathways

426

showed significantly negative correlation with NKA activities with Pearson’s

427

correlation value greater than 0.992 (Table 2).

428

In liver, immune-related KEGG pathways in liver were positively correlated with

429

Na+/K+ concentration, but negatively correlated with antioxidant activities (SI Table

430

S12). In muscle, most immune-related KEGG pathways were negatively correlated 21

431

with K+ concentration. However, the KEGG pathways including complement and

432

coagulation cascades (ko04610) and antigen processing and presentation (ko04612)

433

were the exceptions. (SI Table S13).

434 435

436 437 438 439 440 441

Table 2 Correlations between physiological parameters and immune-related KEGG pathways in gill. KEGG pathway

Na+

K+

NKA

SOD

CAT

MDA

ko04062 Chemokine signaling pathway

-0.727

-0.847

0.695

-0.950

-0.792

0.183

ko04610 Complement and coagulation cascades

-0.087

0.109

-0.995*

0.834

0.970

-0.885

ko04611 Platelet activation

-0.306

-0.486

0.955

-0.983

-0.988*

0.634

ko04612 Antigen processing and presentation

-0.232

-0.038

-0.970

0.744

0.924

-0.943

ko04620 Toll-like receptor signaling pathway

-0.555

-0.382

-0.826

0.466

0.735

-1.000*

ko04621 NOD-like receptor signaling pathway

-0.549

-0.702

0.842

-0.997*

-0.911

0.405

ko04622 RIG-I-like receptor signaling pathway

-0.550

-0.377

-0.828

0.471

0.739

-1.000*

ko04623 Cytosolic DNA-sensing pathway

-0.673

-0.804

0.747

-0.971

-0.836

0.257

ko04624 Toll and Imd signaling pathway

0.008

0.203

-1.000*

0.882

0.989*

-0.836

ko04640 Hematopoietic cell lineage

-0.034

0.161

-0.999*

0.862

0.982

-0.859

ko04650 Natural killer cell mediated cytotoxicity

-0.586

-0.416

-0.803

0.432

0.709

-0.999*

ko04657 IL-17 signaling pathway

0.458

0.276

0.884

-0.562

-0.807

0.995*

ko04658 Th1 and Th2 cell differentiation

0.138

0.329

-0.992*

0.936

1.000*

-0.758

ko04659 Th17 cell differentiation

0.086

0.279

-0.997*

0.916

0.997*

-0.791

ko04660 T cell receptor signaling pathway

-0.022

0.174

-0.999*

0.868

0.984

-0.852

ko04662 B cell receptor signaling pathway

-0.460

-0.278

-0.882

0.560

0.805

-0.996*

ko04664 Fc epsilon RI signaling pathway

-0.809

-0.679

-0.578

0.130

0.454

-0.932

ko04666 Fc gamma R-mediated phagocytosis

-0.856

-0.940

0.528

-0.863

-0.645

-0.028

ko04670 Leukocyte transendothelial migration

-0.333

-0.510

0.947

-0.987

-0.983

0.612

+

+

+

+

+

+

Note: Na , Na concentration; K , K concentration; NKA, Na K -ATPase activity; SOD, SOD activity; CAT, CAT activity; MDA, the MDA contents. * p < 0.05 indicated significant correlation relationship between physiological parameters and KEGG pathways.

4. Discussion In this study, a comparative transcriptomic study along with their physiological

442

and morphological responses of gill, liver, and muscle in yellowfin seabream exposed

443

to different hypo salinities was conducted.

444 445

It is known that activation of NKA was considered to be related to the adaptation of cells in hypoosmotic conditions when fish are acutely exposed to a low salt 22

446

environment (Gui et al., 2016). Measurements of NKA activity, as well as Na+ and K+

447

concentration, could reflect ion transportation and homeostasis of fish body under

448

different salinity environments. In present study, we found marked elevating NKA

449

activity solely in gills in FW group, but no significant differences of Na+ and K+

450

concentration were detected under hypoosmotic stress. The results indicated that

451

enhanced NKA activities may be sufficient to cope with extremely low environmental

452

salinity by preventing the loss of Na+ and K+. Our finding was consistent with the

453

result of Laiz-Carrion et al. (2005), who claimed the raise of NKA activity in gilthead

454

sea bream (Sparus auratus) after acclimation to fresh water. However, these paramters

455

didn’t show any obvious difference in the liver and the muscle, in other words, gill is

456

the most sensitive tissue to osmotic regulation.

457

Oxidative stress is one of the major stress responses when fish is exposed to

458

various environmental stressors (Capela et al., 2016; Zheng et al., 2018a). Monitoring

459

the activity of antioxidant enzymes including SOD and CAT can imply the scavenge

460

capacity of the antioxidant system (Zheng et al., 2019). We examined the SOD, CAT

461

and the MDA content of three tissues in yellowfin seabream and found more obvious

462

responses of these parameters in liver than in gills and muscles in FW group, meaning

463

liver plays a vital role in detoxification. These findings are in agreement

464

with Adeyemi (2014), who claimed oxidative stress were more pronounced in liver

465

compared with that in gill and muscle.

466

In addition, the histological alternation was visualized to present the oxidative

467

stress related morphological changes in tissues (Zheng et al., 2018a, 2018b). Based on 23

468

our observations, more myxocytes, lamellae fusion and hyperplasia in gill were

469

identified in FW group compared with those in BW group, which can be considered

470

as the signs of respiration distress when fish encounters environmental stress

471

(Fernandes et al., 2018; McDonald, 1983). Meanwhile, in FW group, more

472

vacuolations were observed in liver, which is in line with the findings of Fanta et al.,

473

(2003), who had the similar observation in the liver of Corydoras paleatus when they

474

were under exogenous pressure. Whereas, no obvious changes were found in muscles

475

across three treatments (SI Figure S5).

476

Lastly, transcriptome analysis was utilized to understand the immune response and

477

underlying pathways and mechanisms of host to salinity challenges. The results

478

demonstrated that hypoosmotic stress could induce dramatic gene transcription

479

alterations and related KEGG pathway regulations in the gill, liver, and muscle of

480

yellowfin seabream. More osmoregulation and immune related DEGs were identified

481

in gill compared with those in liver and muscle, further supporting that gill is one of

482

the major tissues to maintain osmotic balance by coordinated water and ion transport.

483

The results showed consistency with the finding in many teleost fish (Evans, 2008;

484

Laurent and Dunel, 1980). Several DEGs identified in our study might be the key

485

genes involved in osmoregulation. For example, hsp70 which was found up-regulated

486

17-22 fold in gill of FW group was also detected as a significantly regulated DEG

487

when silver sea bream (Sparus sarba) was subjected to osmotic stress (Deane and

488

Woo, 2004). The other DEG slc6a15 which was identified strikingly up-regulated

489

(92-fold) in liver of FW group was reported as a key gene in liver of spotted sea bass 24

490

for salinity adaptation (Zhang et al., 2017). A set of identified immune-related DEGs

491

in our study such as interleukin and interleukin receptor encoding genes, were also

492

reported as significantly regulated DEGs after salinity challenge in Asian seabass (Xia

493

et al., 2013), spotted sea bass (Zhang et al., 2017), and pipefish (Birrer et al., 2012).

494

And the complement components could assist to link innate and adaptive immunity in

495

fish (Plouffe et al., 2005; Wang and Secombes, 2003; Xia and Yue, 2010). The

496

regulation pattern of DEGs associated with immune response is consistent with those

497

aroused by osmoregulation, indicating that individuals in FW group suffered more

498

stress and had more active immune responses than those in LW group under

499

hypoosmotic stress. Plus, more profound fold changes of DEGs were found in FW

500

group than those in LW group.

501

In addition, the significantly stress-related KEGG pathways were analyzed. The

502

results showed majority of immune-related KEGG pathways were enriched in gill in

503

FW group, followed by liver and then muscle. Again, it demonstrated that yellowfin

504

seabream in in fresh water suffered more stress and had more active immune

505

responses than those in low-saline water. Overall, T cell-mediated immunity was one

506

of the most significant responses when yellowfin seabream was subjected to

507

hypoosmotic stress, which was consistent with the findings that T-cell activation was a

508

key immunological consequence under osmotic stress in striped catfish based on

509

proteomics analysis (Schmitz et al., 2017). Significantly enriched immune related

510

KEGG pathways work together to tackle the exogenous stressor. For example, in gill,

511

antigen-presenting cells (APCs) could process antigens and present antigen peptide 25

512

major histocompatibility complex (MHC) -molecular complexes (ko04612) to CD4+

513

T cells through T cell receptor signaling pathway (ko04660) (Heath and Carbone,

514

2001; Jensen, 2007), then CD4+ T cells are stimulated and differentiated into distinct

515

T helper (Th) cell subsets with different cytokines secretion. To be specific, Th1 cell

516

supports cellular immunity, Th2 cell is functioned as B cell helper, and Th17 cells

517

work against extracellular stressor (Huang et al., 2012). Interleukin-17 (IL-17), a

518

family of pro-inflammatory cytokine, could provoke both acute and chronic

519

inflammatory responses including recruiting of inflammatory cells, producing

520

anti-stressor peptides, and inducing secondary pro-inflammatory mediators by

521

endothelial and epithelial tissues (ko04657). Th17 cell differentiation pathway

522

(ko04659), cooperating with inter-regulating Th1 and Th2 cell differentiation

523

(ko04658) could mediate defensive mechanisms to various stressor (Mosrnann and

524

Coffman, 1989).

525

Competing interests

526 527 528

The authors have no competing interests to declare. CRediT author statement Jianguo Lu: Conceptualization, Methodology, Resources, Writing - Review &

529

Editing, Supervision, Project Administration, Funding Acquisition. Genmei Lin:

530

Methodology, Validation, Formal Analysis, Investigation, Data Curation, Writing -

531

Original Draft, Writing - Review & Editing, Visualization. Min Zheng: Formal

532

Analysis, Investigation, Data Curation, Writing - Review & Editing. Shizhu Li:

533

Methodology. Dong Gao: Methodology. Wenyu Fang: Validation. Jing Huang: 26

534

Investigation. Jingui Xie: Investigation. Jingxiong Liu: Investigation. Yijing Liu:

535

Investigation. Zhaohong Li: Investigation.

536

Acknowledgements

537 538

This work was supported by Guangzhou Science and Technology Project [No. 201803020017] and National Natural Science Foundation of China [No. 31902427].

539

27

540

References

541

Abou-Seedo FS, Dadzie S, Al-Kanaan KA. Sexuality, sex change and maturation patterns in the yellowfin

542

seabream, Acanthopagrus latus (Teleostei: Sparidae) (Houttuyn, 1782). Journal of Applied Ichthyology, 2003,

543

19(2): 65-73. Doi: 10.1046/j.1439-0426.2003.00355.x.

544

Adeyemi JA. Oxidative stress and antioxidant enzymes activities in the African catfish, Clarias gariepinus,

545

experimentally challenged with Escherichia coli and Vibrio fischeri. Fish Physiology and Biochemistry, 2014,

546

40: 347-354. Doi: 10.1007/s10695-013-9847-x.

547

Acar Ü, Parrino V, Kesbiç OS, Lo Paro G, Saoca C, Abbate F, Yılmaz S, Fazio F. Effects of different levels of

548

pomegranate seed oil on some blood parameters and disease resistance against Yersinia Ruckeri in Rainbow

549

Trout. Frontiers in Physiology, 2018, 9: 596. Doi: 10.3389/fphys.2018.00596.

550 551 552 553 554 555 556 557

Birrer SC, Reusch TBH, Roth O. Salinity change impairs pipefish immune defence. Fish & Shellfish Immunology, 2012, 33(6): 1238-1248. Doi: 10.1016/j.fsi.2012.08.028. Boeuf G and Payan P. How should salinity influence fish growth?. Comparative Biochemistry and Physiology Part C: Toxicology & Pharmacology, 2001, 130(4): 411-423. Doi: 10.1016/S1532-0456(01)00268-X. Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics, 2014, 30(15): 2114-2120. Doi: 10.1093/bioinformatics/btu170. Bowden TJ. Modulation of the immune system of fish by their environment. Fish & Shellfish Immunology, 2008, 25(4): 373-383. Doi: 10.1016/j.fsi.2008.03.017.

558

Capela R, Raimundo J, Santos MM, Caetano M, Micaelo C, Vale C, Guimarães L, Reis-Henriques MA. The use of

559

biomarkers as integrative tools for transitional water bodies monitoring in the Water Framework Directive

560

context - A holistic approach in Minho river transitional waters. Science of The Total Environment, 2016, 539:

561

85-96. Doi: 10.1016/j.scitotenv.2015.08.113.

562

Cappello T, Maisano M, Giannetto A, Parrino V, Mauceri A, Fasulo S. Neurotoxicological effects on marine

563

mussel Mytilus galloprovincialis caged at petrochemical contaminated areas (eastern Sicily, Italy): 1H NMR

564

and immunohistochemical assays. Comparative Biochemistry and Physiology Part C: Toxicology &

565

Pharmacology, 2015, 169: 7-15. Doi: 10.1016/j.cbpc.2014.12.006.

566

Cuesta A, Laiz-Carrión R, del Río MPM, Meseguer J, Mancera JM, Esteban MA. Salinity influences the humoral

567

immune parameters of gilthead seabream (Sparus aurata L.). Fish & Shellfish Immunology, 2005,18(3):

568

255-261. Doi: 10.1016/j.fsi.2004.07.009. 28

569

Deane EE and Woo NYS. Differential gene expression associated with euryhalinity in sea bream (Sparus sarba).

570

American Journal of Physiology Regulatory, Integrative and Comparative Physiology, 2004, 287(5):

571

R1054-R1063. Doi: 10.1152/ajpregu.00347.2004.

572

Evans DH. Teleost fish osmoregulation: what have we learned since August Krogh, Homer Smith, and Ancel Keys.

573

American Journal of Physiology Regulatory, Integrative and Comparative Physiology, 2008, 295(2):

574

R704-R713. Doi: 10.1152/ajpregu.90337.2008.

575

Fanta E, Rios FS, Romao S, Vianna ACC, Freiberger S. Histopathology of the fish Corydoras paleatus

576

contaminated with sublethal levels of organophosphorus in water and food. Ecotoxicology and Environmental

577

Safety, 2003, 54: 119-130. Doi: 10.1016/S0147-6513(02)00044-1.

578

Farshadian R, Salati AP, Keyvanshokooh S, Pasha-Zanoosi H. Physiological responses of Yellowfin seabream

579

(Acanthopagrus latus) to acute salinity challenge. Marine and Freshwater Behaviour and Physiology, 2018,

580

51(5): 313-325. Doi: 10.1080/10236244.2019.1573638.

581

Fernandes AL, Nascimento JP, Santos AP, Furtado CA, Romano LA, da Rosa CE, Monserrat JM, Ventura-Lima J.

582

Assessment of the effects of graphene exposure in Danio rerio: A molecular, biochemical and histological

583

approach to investigating mechanisms of toxicity. Chemosphere, 2018, 458-466. Doi:

584

10.1016/j.chemosphere.2018.06.183.

585

Fazio F, Saoca C, Costa G, Zumbo A, Piccione G, Parrino V. Flow cytometry and automatic blood cell analysis in

586

striped bass Morone saxatilis (Walbaum, 1792): A new hematological approach. Aquaculture, 2019, 513:

587

734398. Doi: 10.1016/j.aquaculture.2019.734398.

588

Grabherr MG, Haas BJ, Yassour M, Levin JZ, Thompson DA, Amit I, Adiconis X, Fan L, Raychowdhury R, Zeng

589

Q, Chen Z, Mauceli E, Hacohen N, Gnirke A, Rhind N, di Palma F, Birren BW, Nusbaum C, Lindblad-Toh K,

590

Friedman N, Regev A. Full-length transcriptome assembly from RNA-seq data without a reference genome.

591

Nature Biotechnology, 2011, 29(7): 644-652. Doi: 10.1038/nbt.1883.

592

Gui L, Zhang P, Liang X, Su M, Wu D, Zhang J. Adaptive responses to osmotic stress in kidney-derived cell lines

593

from Scatophagus argus, a euryhaline fish. Gene, 2016, 583: 134-140. Doi: 10.1016/j.gene.2016.02.026.

594

Heath WR and Carbone FR. Cross-presentation in viral immunity and self-tolerance. Nature Reviews Immunology,

595 596 597

2001, 1: 126-134. Doi: 10.1038/35100512. Hu YH, Liu YJ, Tian LX, Yang HJ, Liang GY, Gao W. Optimal dietary carbohydrate to lipid ratio for juvenile Yellowfin Seabream (Sparus latus). Aquaculture Nutrition, 2007, 13(4): 291-297. Doi:

29

598 599 600 601 602 603

10.1111/j.1365-2095.2007.00476.x. Huang G, Wang Y, Chi H. Regulation of TH17 cell differentiation by innate immune signals. Cellular and Molecular Immunology, 2012, 9: 287-295. Doi: 10.1038/cmi.2012.10. Jensen PE. Recent advances in antigen processing and presentation. Nature Immunology, 2007, 8(10): 1041-1048. Doi: 10.1038/ni1516. Jiang I-F, Bharath Kumar V, Lee D-N, Weng C-F. Acute osmotic stress affects Tilapia (Oreochromis mossambicus)

604

innate immune responses. Fish & Shellfish Immunology, 2008, 25(6): 841-846. Doi: 10.1016/j.fsi.2008.09.

605

006.

606

Klopfenstein DV, Zhang L, Pedersen BS, Ramírez F, Vesztrocy AW, Naldi A, Mungall CJ, Yunes JM, Botvinnik O,

607

Weigel M, Dampier W, Dessimoz C, Flick P, Tang H. GOATOOLS: A Python library for Gene Ontology

608

analyses. Scientific Reports, 2018, 8: 10872. Doi: 10.1038/s41598-018-28948-z.

609

Laiz-Carrion R, Gurreiro PM, Fuentes J, Canario AVM, Martin Del Rio MP, Mancera JM. Branchial

610

osmoregulatory response to salinity in the gilthead sea bream, Sparus auratus. Jorunal of Experimental

611

Zoology, 2005, 303A(7): 563-576. Doi: 10.1002/jez.a.183.

612

Laurent P and Dunel S. Morphology of gill epithelia in fish. American Journal of Physiology: Regulatory,

613

Integrative and Comparative Physiology, 1980, 238: R147-R159. Doi: 10.1152/ajpregu.1980.238.3.R147.

614

Li B and Dewey CN. RSEM: accurate transcript quantification from RNA-Seq data with or without a reference

615 616 617 618 619 620

genome. BMC Bioinformatics, 2011,12:323. Doi: 10.1186/1471-2105-12-323. Li W and Godzik A. Cd-hit: a fast program for clustering and comparing large sets of protein or nucleotide sequences. Bioinformatics, 2006, 22(13): 1658-1659. Doi: 10.1093/bioinformatics/btl158. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biology, 2014, 15(12): 550. Doi: 10.1186/s13059-014-0550-8. Mao X, Cai T, Olyarchuk JG, Wei L. Automated genome annotation and pathway identification using the KEGG

621

Orthology (KO) as a controlled vocabulary. Bioinformatics, 2005, 21(19): 3787-3793, Doi:

622

10.1093/bioinformatics/bti430.

623

Marc AM, Quentel C, Severe A, Le Bail PY, Boeuf G. Changes in some endocrinological and non-specific

624

immunological parameters during seawater exposure in the brown trout. Journal of Fish Biology, 1995, 46(6):

625

1065-1081. 10.1111/j.1095-8649.1995.tb01410.x.

626

McDonald DG. The effects of H+ upon the gills of freshwater fish. Canadian Journal of Zoology, 1983, 61(4):

30

627 628

691-703. Doi: 10.1139/z83-093. Mosrnann TR and Coffman RL. TH1 and TH2 cells: different patterns of lymphokine secretion lead to different

629

functional properties. Annual Review of Immunology, 1989, 7: 145-173. Doi: 10.1146/annurev.iy.07.040189.

630

001045.

631

Movahedinia AA, Savari A, Morovvati H, Koochanin P, Marammazi JG, Nafisi M. The effects of changes in

632

salinity on gill mitochondria-rich cells of Juvenile Yellowfin Seabream, Acanthopagrus latus. Journal of

633

Biological Sciences, 2009, 9(7): 710-720. Doi: 10.3923/jbs.2009.710.720.

634 635

Nakabo T. Fishes of Japan: with pictorial keys to the species, English edition. Tokai University Press, 2002, Kanagawa.

636

Parrino V, Cappello T, Costa G, Cannavà C, Sanfilippo M, Fazio F, Fasulo S. Comparative study of haematology

637

of two teleost fish (Mugil cephalus and Carassius auratus) from different environments and feeding habits.

638

The European Zoological Journal, 2018, 85(1): 193-199. Doi: 10.1080/24750263.2018.1460694.

639

Parrino V, Costa G, Cannavà C, Fazio E, Bonsignore M, Concetta S, Piccione G, Fazio F. Flow cytometry and

640

micro-Raman spectroscopy: Identification of hemocyte populations in the mussel Mytilus galloprovincialis

641

(Bivalvia: Mytilidae) from Faro Lake and Tyrrhenian Sea (Sicily, Italy). Fish and Shellfish Immunology,

642

2019a, 87: 1-8. Doi: 10.1016/j.fsi.2018.12.067.

643

Parrino V, Kesbiç OS, Acar Ü, Fazio F. Hot pepper (Capsicum sp.) oil and its effects on growth performance and

644

blood parameters in rainbow trout (Oncorhynchus mykiss). Natural Product Research, 2019b: 1-5. Doi:

645

10.1080/14786419.2018.1550769.

646 647 648 649

Plouffe DA, Hanington PC, Walsh JG, Wilson EC, Belosevic M. Comparison of select innate immune mechanisms of fish and mammals. Xenotransplantation, 2005, 12(4): 266-277. Doi: 10.1111/j.1399-3089.2005.00227.x. Rizopoulos D. ltm: An R package for latent variable modelling and item response theory analyses. Journal of Statistical Software, 2006, 17(5): 1-25. http://www.jstatsoft.org/v17/i05/

650

Schmitz M, Douxfils J, Mandiki SNM, Morana C, Baekelandt S, Kestemont P. Chronic hyperosmotic stress

651

interferes with immune homeostasis in striped catfish (Pangasianodon hypophthalmus, S.) and leads to

652

excessive inflammatory response during bacterial infection. Fish & Shellfish Immunology 2016, 55: 550-558.

653

Doi: 10.1016/j.fsi.2016.06.031.

654 655

Schmitz M, Ziv T, Admon A, Baekelandt S, Mandiki SNM, L'Hoir M, Kestemont P. Salinity stress, enhancing basal and induced immune responses in striped catfish Pangasianodon hypophthalmus (Sauvage). Journal of

31

656 657

Proteomics, 2017, 167(7): 21-24. Doi: 10.1016/j.jprot.2017.08.005. Simão FA, Waterhouse RM, Ioannidis P, Kriventseva EV, Zdobnov EM. BUSCO: assessing genome assembly and

658

annotation completeness with single-copy orthologs. Bioinformatics, 2015, 31(19): 3210-3212. Doi:

659

10.1093/bioinformatics/btv351.

660 661 662 663 664

Smith-Unna R, Boursnell C, Patro R, Hibberd JM, Kelly S. TransRate: reference-free quality assessment of de novo transcriptome assemblies. Genome Research, 2016, 26: 1134-1144. Doi: 10.1101/gr.196469.115. Storey JD, Bass AJ, Dabney A, Robinson D. qvalue: Q-value estimation for false discovery rate control. R package version, 2019, 2.18.0. http://github.com/jdstorey/qvalue Su M, Zhou J, Duan Z, Zhang J. Transcriptional analysis of renal dopamine-mediated Na+ homeostasis response to

665

environmental salinity stress in Scatophagus argus. BMC Genomics, 2019, 20: 418. Doi:

666

10.1186/s12864-019-5795-x.

667

Taylor JF, Needham MP, North BP, Morgan A, Thompson K, Migaud H. The influence of ploidy on saltwater

668

adaptation, acute stress response and immune function following seawater transfer in non-smolting rainbow

669

trout. General and Comparative Endocrinology, 2007, 152(2-3): 314-325. Doi: 10.1016/j.ygcen.2007.02.029.

670 671 672

Tort L. Stress and immune modulation in fish. Developmental & Comparative Immunology, 2011, 35(12): 1366-1375. Doi: doi.org/10.1016/j.dci.2011.07.002. Vargas-Chacoff L, Saavedra E, Oyarzún R, Martínez-Montaño E, Pontigo JP, Yáñez A, Ruiz-Jarabo I, Mancera JM,

673

Ortiz E, Bertrán C. Effects on the metabolism, growth, digestive capacity and osmoregulation of juvenile of

674

Sub-Antarctic Notothenioid fish Eleginops maclovinus acclimated at different salinities. Fish Physiology and

675

Biochemistry, 2015, 41(6): 1369-1381. Doi: 10.1007/s10695-015-0092-3.

676

Wang T and Secombes CJ. Complete sequencing and expression of three complement components, C1r, C4 and C1

677

inhibitor, of the classical activation pathway of the complement system in rainbow trout Oncorhynchus mykiss.

678

Immunogenetics, 2003, 55(9): 615-628. Doi: 10.1007/s00251-003-0622-5.

679

Wu XY, Liu YJ, Tian LX, Mai KS, Yang HJ, Liang GY. Effects of raw corn starch levels on growth, feed

680

utilization, plasma chemical indices and enzyme activities in juvenile yellowfin seabream Sparus latus

681

Houttuyn. Aquaculture Research, 2007a, 38(12): 1330-1338. Doi: 10.1111/j.1365-2109.2007.01812.x.

682

Wu XY, Liu YJ, Tian LX, Mai KS, Yang HJ. Utilization of different raw and pre-gelatinized starch sources by

683

juvenile yellowfin seabream Sparus latus. Aquaculture Nutrition, 2007b, 13(5): 389-396. Doi:

684

10.1111/j.1365-2095.2007.00489.x.

32

685

Xia JH, Huang JH, Gong JB, Jiang SG. Significant population genetic structure of yellowfin seabream

686

Acanthopagrus latus in China. Journal of Fish biology, 2008, 73(8): 1979-1992. Doi:

687

10.1111/j.1095-8649.2008.02010.x.

688

Xia JH, Liu P, Liu F, Lin G, Sun F, Tu R, Yue GH. Analysis of stress-responsive transcriptome in the intestine of

689

Asian seabass (Lates calcarifer) using RNA-Seq. DNA Research, 2013, 20(5): 449-460. Doi:

690

10.1093/dnares/dst022.

691 692 693 694 695

Xia JH, Xia KF, Jiang SG. Characterization of 11 polymorphic microsatellite loci in the yellowfin seabream Acanthopagrus latus. Molecular Ecology Notes, 2006, 6: 484-486. Doi: 10.1111/j.1471-8286.2006.01282.x. Xia JH and Yue GH. Identification and analysis of immune-related transcriptome in Asian seabass Lates calcarifer. BMC Genomics, 2010, 11: 356. Doi: 10.1186/1471-2164-11-356. Zhang X, Wen H, Wang H, Ren Y, Zhao J, Li Y. RNA-Seq analysis of salinity stress-responsive transcriptome in

696

the liver of spotted sea bass (Lateolabrax maculatus). PLoS ONE, 2017, 12(3): e0173238. Doi:

697

10.1371/journal.pone.0173238.

698

Zheng M, Lu J, Zhao D. Effects of starch-coating of magnetite nanoparticles on cellular uptake, toxicity and gene

699

expression profiles in adult zebrafish. Science of The Total Environment, 2018a, 622-623: 930-941. Doi:

700

10.1016/j.scitotenv.2017.12.018.

701

Zheng M, Lu J, Zhao D. Toxicity and transcriptome sequencing (RNA-seq) analyses of adult zebrafish in response

702

to exposure carboxymethyl cellulose stabilized iron sulfide nanoparticles. Scientific reports, 2018b, 8: 8083.

703

Doi: 10.1038/s41598-018-26499-x.

704

Zheng M, Lu J, Lin G, Su H, Sun J, Luan T. Dysbiosis of gut microbiota by dietary exposure of three

705

graphene-family materials in zebrafish (Danio rerio). Environmental Pollution, 2019, 254: 112969. Doi:

706

10.1016/j.envpol.2019.112969.

33

Highlights: ·First tissue-specific RNA-seq analysis on hypoosmotic stress in A. latus. ·Fish in fresh water suffered more stress than those in low-saline water. ·Gill harbored the most osmoregulation and immune-related DEGs. ·T cell-mediated immunity pathways were essential to tackle hypoosmotic stress.