Transcriptome profiling of grass carp (Ctenopharyngodon idellus) infected with Aeromonas hydrophila

Transcriptome profiling of grass carp (Ctenopharyngodon idellus) infected with Aeromonas hydrophila

Accepted Manuscript Transcriptome profiling of Grass carp (Ctenopharyngodon idellus) infected with Aeromonas hydrophila Ying Yang, Hui Yu, Hua Li, Anl...

738KB Sizes 0 Downloads 45 Views

Accepted Manuscript Transcriptome profiling of Grass carp (Ctenopharyngodon idellus) infected with Aeromonas hydrophila Ying Yang, Hui Yu, Hua Li, Anli Wang PII:

S1050-4648(16)30079-1

DOI:

10.1016/j.fsi.2016.02.035

Reference:

YFSIM 3847

To appear in:

Fish and Shellfish Immunology

Received Date: 31 December 2015 Revised Date:

25 February 2016

Accepted Date: 29 February 2016

Please cite this article as: Yang Y, Yu H, Li H, Wang A, Transcriptome profiling of Grass carp (Ctenopharyngodon idellus) infected with Aeromonas hydrophila, Fish and Shellfish Immunology (2016), doi: 10.1016/j.fsi.2016.02.035. This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. 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.

ACCEPTED MANUSCRIPT 1

Transcriptome profiling of Grass carp (Ctenopharyngodon idellus) infected with

2

Aeromonas hydrophila

3

Ying Yanga, b, Hui Yub**, Hua Lib, Anli Wanga*

4 a

College of Life Sciences, South China Normal University, Guangzhou, Guangdong 510631, China

6 7

b

College of Life Science, Foshan University, Foshan, Guangdong 528231, China

8

SC

9 10

M AN U

11

*Correspondence to:

13

Pro. Anli Wang

14

College of Life Science

15

South China Normal University

16

No. 55 West Zhongshan Avenue

17

Tianhe, Guangzhou, Guangdong, P.R.China

18

Tel/fax +86 2085210141.

19

E-mail: [email protected]

TE D

12

20

**Co- Correspondence to:

22

Pro. Hui Yu

23

College of Life Science

24

Foshan University

25

No.1 Xianhu Road

26

Nanhai, Foshan, Guangdong, P.R.China

27

Tel/fax +86 757 86678587.

28

E-mail: [email protected]

30 31

AC C

EP

21

29

RI PT

5

ACCEPTED MANUSCRIPT Abstract: Aeromonas hydrophila is the causative pathogen of intestinal hemorrhage

33

which has caused great economic loss in grass carp aquaculture. In order to

34

understand the immunological response of grass carp to infection by A. hydrophila,

35

the transcriptomic profiles of the spleens from infected and non-infected grass carp

36

groups were obtained using HiSeq™ 2500 (Illumina). An average of 63 million clean

37

reads per library was obtained, and approximately 80% of these genes were

38

successfully mapped to the reference genome. A total of 1591 up-regulated and 530

39

down-regulated genes were identified. Eight immune-related categories involving 105

40

differently expressed genes were scrutinized. 16 of the differently expressed genes

41

involving immune response were further validated by qRT-PCR. Our results provide

42

valuable information for further analysis of the mechanisms of grass carp defense

43

against A. hydrophila invasion.

44

Key words: Ctenopharyngodon idellus; Aeromonas hydrophila; Transcriptome;

45

Pathway analysis; Immune response

46

1. Introduction

47

Grass carp (Ctenopharygodon idella) is an economically important species and its

48

global production is more than 4.5 million tones per year, which making it the most

49

highly consumed freshwater fish worldwide[1]. However, infectious intestinal

50

hemorrhage caused by Aeromonas hydrophila has been severe for years[2], resulting

51

in great economic loss and threatening the development of grass carp aquaculture. A.

52

hydrophila is facultative anaerobic, motile, and Gram-negative rods in the family

53

Aeromonadaceae[3]. It is often found in association with hemorrhagic septicemia in

54

cold-blooded animals including fish, reptiles and amphibians. However, this organism

55

has also attracted attention as an emerging human pathogen, and has been considered

56

to have a significant impact on public health[4]. Therefore, identification of host

57

factors in response to A. hydrophila infection has a great significance for disease

58

prevention and control in grass carp culture.

AC C

EP

TE D

M AN U

SC

RI PT

32

59

Transcriptome analysis is a powerful tool for leading to a better understanding of

60

the underlying pathways and mechanisms controlling cell fate, development, and

61

disease progression[5]. Over the years, several technologies have been developed to

ACCEPTED MANUSCRIPT survey transcriptomes in a high-throughput manner. High-throughput sequencing

63

(HTS) technologies permit genome-wide transcriptomic analysis at a higher resolution,

64

and these technologies have been widely used to study pathogenic processes during

65

bacterial infections[6]. Hegedus et al. firstly used Solexa/Illumina’s DGE system to

66

study zebrafish transcriptome after Mycobacterium marinum infection, it showed the

67

feasibility and superiority of high-throughput sequencing technology in the field of

68

fish immune research[7]. Then RNA-Seq have been applied to immune-related gene

69

and signaling pathway analysis of several fish species such as Mandarin fish, Large

70

yellow croaker and Lates calcarijer[8-10].

SC

RI PT

62

The results of sequencing of the transcriptome of juvenile grass carp after A.

72

hydrophila infection and non-infection were presented in this study. High-quality

73

cDNA sequences were obtained by using Illumina RNA-Seq method and the

74

sequencing reads were mapped to the reference genome database of grass carp[11].

75

Additionally, a great number of immune related genes that were differently expressed

76

upon A. hydrophila infection were obtained and functionally annotated, and the gene

77

expression patterns of some of these genes were verified by qRT-PCR. These results

78

provide a valuable resource for further research into the mechanism of anti-infection

79

immunity of grass carp.

80

2. Materials and methods

81

2.1 Fish and bacteria

82

Healthy juvenile grass carps of 20g±2g body weight were kindly provided by the

83

Hold-one aquatic breeding center in China and reared under pathogen-free conditions.

84

Fish were maintained in aerated tap water at 28 °C in aquaria with Eheim biofilters

85

until use. The bacterial strain A. hydrophila used for the experiments was isolated

86

from diseased grass carp, and has been deposited in the China Center for Type Culture

87

Collection under preservation number CCAM05068. The bacteria was cultured in

88

Luria-Bertain (LB) at 28 °C for 24 h with constant shaking (150 rpm), then harvested

89

by centrifugation at 6000 rpm for 5 min, washed once by PBS and centrifuged again.

90

Bacterial pellets were then re-suspended in saline adjusted the concentration to 5×107

91

colony forming units (CFU) ml-1.

AC C

EP

TE D

M AN U

71

ACCEPTED MANUSCRIPT 2.2 Challenge experiments and RNA preparation

93

Healthy fish were randomly divided into two groups (30 fish per group). For bacterial

94

infection group, fish were injected intraperitoneally (i.p.) with 0.1ml A. hydrophila

95

suspension above-mentioned. The fish treated with 0.1ml 0.65% physiological saline,

96

were used as the control. After 6h of injection, the spleen were collected and

97

immediately stored in liquid nitrogen until RNA extraction.

RI PT

92

Total RNA was extracted from each sample using TRIzol Reagent (Invitrogen,

99

USA) according to the manufacture's protocol. The quality and quantity of the RNA

100

samples were examined by use of the Agilent 2100 Bioanalyzer (Agilent Technologies)

101

and the integrity was assessed by electrophoresis on 1% agarose gel. For each group,

102

equal amount of RNA from the nine fish individuals per line were pooled to provide

103

templates for RNA-Seq library construction.

104

2.3 Library construction and Illumina sequencing

105

After DNase I treatment, mRNA was purified using oligo (dT)25 magnetic beads

106

(Dynabeadsoligo (dT)25, Invitrogen) and subsequently interrupted to short

107

fragments of 200-250 nucleotides using RNA fragmentation reagent (Ambion, USA).

108

Then sequencing libraries were generated using NEBNext® Ultra™ RNA Library

109

Prep Kit for Illumina® (NEB, USA) following manufacturer’s recommendations. The

110

cDNAs were checked by Agilent 2100 Bioanalyzer (Agilent, USA) and ABI

111

StepOnePlus Real-Time PCR System (ABI, USA). The mixed DNA libraries were

112

diluted to 4-5 pM for sequencing by Illumina Hiseq 2500™ platform.

113

2.4 Sequence annotation

114

Raw reads were first cleaned by removing adaptor sequences, low quality sequences

115

(Sanger base quality < 20) and reads with unknown nucleotides larger than 10%.

116

Clean reads were mapped to the Ctenopharyngodon idellus reference genome

117

(http://www.ncgr.ac.cn/grasscarp/) using HISAT (version 0.1.6) [12]. The HTSeq

118

(version 0.6.1) [13] was utilized to calculate the number of aligned reads per exon

119

through annotation of the grass carp genome. Subsequently, the transcripts were

120

subjected to BLASTX similarity searching against NCBI non-redundant protein

121

database (NR) with an E-value threshold of 10-5, and the unigenes were identified.

AC C

EP

TE D

M AN U

SC

98

ACCEPTED MANUSCRIPT 2.5 Differential expression analysis

123

Expression levels were measured in reads per kilo base of exon per million mapped

124

reads (RPKM) method. The distribution of gene expression was analyzed. Statistical

125

comparison between two groups was conducted using DEGseq (version 1.18.0) [14].

126

To assess the significance of differential gene expression, the threshold of FDR was

127

set at ≤ 0.01 and the absolute value of log2 ratios (fold change between non-infected /

128

infected samples) at ≥1. Differential expression genes (DEGs) were further annotated

129

by Gene ontology (GO) functional enrichment and Kyoto Encyclopedia of Genes and

130

Genomes (KEGG) pathway analysis. GO annotation was acquired using the GOseq

131

(version 1.16.2) and WEGO [15, 16]. P-values generated from the enrichment

132

analysis were subjected to multiple hypothesis testing, with p-values ≤ 0.05

133

considered significant. Pathway analysis was performed using the Kyoto

134

Encyclopedia

135

(http://www.genome.jp/kegg)[17]. After multiple test correction, pathways with

136

Q-values ≤ 0.05 were considered to be significantly enriched in DEGs.

137

2.6 Confirmation using quantitative real-time RT-PCR (QPCR)

138

To examine the reliability of the RNA-Seq results, a selected subset of differently

139

expressed genes involved in immune responses were selected for validation using

140

quantitative real-time RT-PCR (qRT-PCR). 16 genes differently expressed between

141

the infected and non-infected groups, including Chemokine C motif receptor 1 (XCR

142

1), interleukin-1 beta (IL-1β), tumor necrosis factorα (TNF-α), C-type lysozyme

143

(LysC), interferon regulatory factor 4 (IRF4), B cell receptor CD22, coagulation

144

factor 7 (F7), complement component 1 (C1), complement component 3 (C3),

145

complement component 7 (C7), toll-like receptor 22 (TLR 22), MHC class I antigen

146

(MHC-I Ag), MHC class II antigen (MHC-II Ag), macrophage expressed gene 1

147

(MPEG 1), heat shock 70 kDa protein (HSP 70), and α-2-macroglobulin (α2MG)

148

were selected for qRT-PCR assay. The housekeeping gene β-actin was used as the

149

reference gene. Suitable primers were designed using Primer Express 3.0(Table.1) and

150

synthesized by Sangon Bothech (Shanghai) Co., Ltd. QRT-PCR with SYBR Green

151

dye (TaKaRa, Dalian, China) was performed on an ABI PRISM 7500 Fast Real-Time

Genes

and

Genomes

(KEGG)

database

AC C

EP

TE D

of

M AN U

SC

RI PT

122

ACCEPTED MANUSCRIPT PCR System according to the manufacturer’s protocol. Each treatment group was

153

created by combining equal amounts of RNA from three replicate pools (three

154

individual fish per pool). All reactions were performed in triplicates. The qRT-PCR

155

conditions were as follows: 30s at 95, followed 40 cycles of 5s at 95, 30s at 58,

156

34s at 72. Dissociation curve analysis was performed to determine the target

157

specificity. The relative expression ratio of the target genes versus β-actin gene was

158

calculated using 2-△△CT method and all data were given in terms of relative mRNA

159

expression.

RI PT

152

SC

Insert Table 1 here.

160

3. Results

162

3.1 Transcriptome sequence and reads mapping

163

Infected and non-infected groups were analyzed by RNA-Seq respectively. RNA

164

sequencing resulted in about 65.5 million (infected) and 63.4 million (non-infected)

165

raw reads of 125 bp (Table 2). After filtration, 64.6 million of infected groups and

166

62.6 million of non-infected groups clean reads were obtained, and these sequences

167

were carried forward for additional analysis. Average 80% of the clean reads were

168

mapped to the reference genome, 24507 genes of infected group and 24383 genes of

169

non-infected group were detected. These sequences were compared with the NCBI

170

non-redundant (NR) protein for functional annotation. Raw sequencing reads data has

171

been submitted to Sequence Read Archive in NCBI, the SRA accession numbers are

172

SRR3045340 and SRR3045341.

TE D

EP

AC C

173

M AN U

161

Insert Table 2 here.

174

3.2 Differential expression analysis

175

RNA-Seq reads mapped to the grass carp reference genome were aligned and their

176

relative abundances were estimated using HTseq. According to the statistical analysis

177

of unigene data with DEGseq, the genes significantly differentially expressed between

178

infected and non-infected groups were identified. The results showed that 2121 genes

179

were significantly differently expressed, including 1591 up-regulated genes and 530

180

down-regulated genes. Differential expression genes and no differential expression

181

genes distribution trends were presented clearly (Fig. 1).

ACCEPTED MANUSCRIPT Insert Fig. 1 here.

182

3.3 GO and KEGG enrichment analysis of differently expressed genes

184

The results of the GO enrichment analysis of differential expression genes were

185

classified into three categories: biological process (121 subclasses), cellular

186

component (31 subclasses) and molecular function (129 subclasses). The significant

187

GO terms in these three categories and DEGs of them were showed in Fig. 2. The

188

pathway enrichment analysis helped to further understand genes’ biological functions.

189

In this study, 1631 of 2121 DEGs were annotated to 251 signaling pathways in KEGG.

190

The top30 most enriched pathway terms and DEGs of them were showed in Fig. 3. Insert Fig. 2 here.

192

Insert Fig. 3 here.

M AN U

191

SC

RI PT

183

3.4 Analysis of differently expressed genes related to immune responses

194

In this report, we focused on the genes related to the immune responses. Based on

195

gene clustering and pathway analysis, genes involved in both innate and adaptive

196

immunity had been identified (Table 3), including 85 up-regulated and 20

197

down-regulated genes. They were further grouped into 8 sub-classes as follows:

198

phagocytosis, complement system, cytokine, antigen processing and presentation,

199

pattern recognition receptor, cell adhesion molecules, apoptosis, antioxidant enzymes.

200

Insert Table 3 here.

EP

TE D

193

3.5 Analysis of differently expressed genes by qRT-PCR

202

Quantitative real-time RT-PCR was performed on 16 immune genes for validating the

203

differently expressed genes identified by RNA-Seq. For the results, the analysis of

204

qRT -PCR demonstrated a single product for all tested genes. Fold changes of qRT

205

-PCR were compared with the RNA-Seq expression profiles. As shown in Fig. 4, all

206

of qRT -PCR results were significantly correlated with the RNA-Seq results and they

207

showed the identical up-regulated and down-regulated patterns of these genes in both

208

assays.

AC C

201

Insert Fig. 4 here.

209 210

4. Discussion

211

In the present study, transcriptome profiling of grass carp was undertaken to gain

ACCEPTED MANUSCRIPT insights into the molecular mechanisms underpinning host responses to A. hydrophila.

213

There were a total of 2121 genes expressed differently between infected and

214

non-infected group and these differently expressed genes were assigned multiple

215

potential functions in molecular function and biological process. Here, we mainly

216

focused on the immune-related genes of fish activated in the early stage of

217

bacteria invasion. Taking the advantages of GO annotation and KEGG pathway

218

classification, we categorized these immune related genes into eight categories and

219

the most important parts of them were discussed.

220

4.1 Antigen processing and presentation

SC

RI PT

212

Apparently, innate immunity is the major approach to defense against invading

222

virulent bacteria like A. hydrophila. But from the path-ways enriched in this present

223

study, cellular immunity belonging to acquired immunity also plays a role in the

224

disease. Based on analysis of gene differential expression, T-lymphoma invasion and

225

metastasis-inducing protein 1, T cell receptor beta chain and MHC class I antigen

226

were up-regulated significantly. CD2 family receptor and other CDs involved in

227

lymphocyte recognition and signal transduction were found in DEGs. In addition, the

228

discovery of novel immune-type receptors (NITRs) indicates that fish leukocytes

229

express various types of Ig superfamily (IgSF) receptors for contact with other

230

immune cells[18]. It is worth noting that MHC I and MHC II were both DEGs, but

231

MHC I was up-regulated in this tissue while MHC II was down-regulated.

232

We noticed that B-cell receptor CD22 was up-regulated and the expression of this

233

membrane receptor is generally regarded as being restricted to B cells[19]. On the

234

contrary, activated cytotoxic T lymphocytes (CTLs) released toxic particles such as

235

perforin to kill host cells infected with pathogens[20].

236

4.2 Phagocytosis

AC C

EP

TE D

M AN U

221

237

Phagocytosis is an important, evolutionarily conserved mechanism that is integral

238

to host defenses against invading microorganisms[21]. In professional phagocytes

239

such as neutrophils and macrophages, phagocytosis is triggered by the recognition of

240

pathogen

241

immunoglobulins and complement components bound to the microbes, which are

associated

molecular

patterns

(PAMPs)

and

opsonins

such

as

ACCEPTED MANUSCRIPT 242

recognized

by

243

differently expressed genes, 23 phagocytosis–related genes were found and all of

244

them were up-regulated. MPEG 1 and neutrophil cytosolic factor were up-regulated

245

indicating both macrophages and neutrophils were stimulated. Calreticulin and

246

calmegin

247

Phagocytosis is initiated by the interaction of receptors on the surface of the

248

phagocyte with ligands on large particles such as bacterial. Subsequently integrins

249

mediate particle internalization[25]. Once internalized, the phagosome that encloses

250

the particle undergoes a series of maturation steps and culminates in phagolysosome

251

fusion[26]. In this study, increasing expression of lysosome-related genes including

252

C-type lysozyme, G-type lysozyme, Cathepsin K, Cathepsin B precursor, Cathepsin L

253

precursor, and Cysteine proteinase promotes the lysis of bacteria[27-30].

254

4.3 Complement system

receptors[22].

enhance the immune competence

of

Among

macrophages[23,

the

24].

M AN U

SC

RI PT

could

phagocytosis-promoting

Complement system is a key component of innate immunity in teleosts and its

256

activation leads to opsonization of pathogens, recruitment of inflammatory and

257

immune competent cells, and the direct killing of pathogens[31, 32]. There are three

258

different ways to activate the complement on the surface of invading pathogens:

259

classical, lectin and alternative pathways[33]. In our results, C1

260

up-regulated, indicating the activation of classical pathway. C3 was down-regulated

261

and C3b, factor B were up-regulated, deducing C3 molecules were decomposed and

262

activate the alternative pathway.

EP

C7, C8β were

The blood coagulation pathway controls the coagulation and fibrinolysis of blood

AC C

263

TE D

255

264

through a series of sequentially activated serine kinases[34]. In our data,

265

mannose-binding protein-associated serine protease, α-2-macroglobulin and several

266

coagulation factors were significantly down-regulated, indicating the inhibition of

267

blood coagulation[35, 36]. This may be related to the bleeding symptoms observed in

268

grass carp caused by A. hydrophila[37].

269

4.4 Cytokines

270

Cytokines are a family of low molecular weight proteins secreted by activated

271

immune-related cells upon induction by various pathogens such as parasitic, bacterial,

ACCEPTED MANUSCRIPT or viral components[38]. Cytokines could be divided into interferons (IFNs),

273

interleukins (ILs), tumor necrosis factors (TNFs), colony stimulating factors, and

274

chemokines[39]. In this report, a total of 20 immune related cytokine genes were

275

found in DEGs, 17 were up-regulated and three were down-regulated. These

276

differently expressed cytokine genes could be classified as interferons and signaling

277

factors (IRF 4, IFN transmembrane protein 1, IRF 10 and IRG-47), interleukin family

278

members and receptors (IL-1β, I IL -10, IL -21 receptor, IL -1 receptor I, IL -27

279

subunit beta precursor and IL 13 receptor precursor), tumor necrosis factors (TNF-α

280

induced protein 6 and TNF-α receptor 2) and chemokines (Chemokine-like receptor 1,

281

XCR1, CXCR3, CXC motif chemokine 9-like and CXC motif chemokine 13

282

precursor). Interleukin, tumor necrosis factors and chemokines play key roles in

283

inflammation and host defense, besides that, activity of IFN and IRF genes expression

284

also proved that interferons involved in bacterial immune responses[40-42]. The

285

complicate coordinated approaches of the cytokine-related genes in the immune

286

responses of grass carp during bacteria infection remains to be addressed in the future.

287

4.5 Others

TE D

M AN U

SC

RI PT

272

288

The results also demonstrate that a number of other differentially expressed

289

immune-related genes enriched in pattern recognition receptor (PRR), cell adhesion

290

molecules

291

previous reports suggested that toll-like receptor 22 (TLR22) exists exclusively in

292

aquatic animals and recognizes double stranded RNA (dsRNA)[43, 44], but this study

293

shows it can also be elevated by bacterial stimulation. The details of the regulation of

294

these immune-related genes need to be elucidated in the future.

295

5. Conclusion

296

In conclusion, RNA-seq analysis is a highly valuable resource to address the issues of

297

disease resistance at the genome and transcriptome levels. This study provided a first

298

survey of host defense gene activities against bacterial challenge based on

299

differentially

300

significantly enriched gene clusters including phagocytosis, complement system,

301

cytokines, antigen processing and presentation, PRR, CAM, apoptosis and antioxidant

apoptosis

and

antioxidant

enzymes.

Remarkably,

AC C

EP

(CAM),

transcriptomic

profiling

in

grass

carp.

Additionally,

ACCEPTED MANUSCRIPT 302

enzymes were identified and discussed. This study will add new information that may

303

help to prevent the A. hydrophila infection.

304

Acknowledgments

306

We would like to thank Yu-dong Du (Project manager of Hold-one aquatic breeding

307

co., LTD) for providing juvenile grass carps. This work was jointly supported by

308

National Spark Program Project of China (2015GA780004), Guangdong Natural

309

Science Foundation (2013B020307017; 2013B020503070) and Foshan Science and

310

Technology Innovation Platform Construction Program (2013AG10005).

SC

RI PT

305

311

References

313

[1] Stephan M, Hobsdawn P. Australian fisheries and aquaculture statistics 2013.

314

Fisheries Research and 2014.

M AN U

312

[2] Liu L, Gong Y-X, Zhu B, Liu G-L, Wang G-X, Ling F. Effect of a new

316

recombinant Aeromonas hydrophila vaccine on the grass carp intestinal

317

microbiota and correlations with immunological responses. Fish Shellfish

318

immun. 2015; 45(1):175-83.

TE D

315

[3] Seshadri R, Joseph SW, Chopra AK, Sha J, Shaw J, Graf J, et al. Genome

320

sequence of Aeromonas hydrophila ATCC 7966T: jack of all trades. J bacteriol.

321

2006; 188(23):8272-82.

323 324

[4] Agger WA, McCormick J, Gurwith MJ. Clinical and microbiological features of Aeromonas

hydrophila-associated

AC C

322

EP

319

diarrhea.

J

Clin

Microbiol.

1985;

21(6):909-13.

325

[5] Mutz KO, Heilkenbrinker A, Lonne M, Walter J-G, Stahl F. Transcriptome

326

analysis using next-generation sequencing. Curr opin biotech. 2013;

327 328 329

24(1):22-30. [6] Soon WW, Hariharan M, Snyder MP. High-throughput sequencing for biology and medicine. Mol syst biol. 2013; 9(1):640.

330

[7] Hegedűs Z, Zakrzewska A, Ágoston VC, Ordas A, Rácz P, Mink M, et al. Deep

331

sequencing of the zebrafish transcriptome response to mycobacterium infection.

ACCEPTED MANUSCRIPT 332

Mol Immunol. 2009; 46(15):2918-30.

333

[8] Khoo CK, Abdul-Murad AM, Kua BC, Mohd-Adnan A. Cryptocaryon irritans

334

infection induces the acute phase response in Lates calcarifer: a transcriptomic

335

perspective. Fish Shellfish immun. 2012; 33(4):788-94. [9] Mu Y, Li M, Ding F, Ding Y, Ao J, Hu S, et al. De novo characterization of the

337

spleen transcriptome of the large yellow croaker (Pseudosciaena crocea) and

338

analysis of the immune relevant genes and pathways involved in the antiviral

339

response. Plos One. 2014; 9(5):e97471.

RI PT

336

[10] Zhou W, Zhang Y, Wen Y, Ji W, Zhou Y, Ji Y, et al. Analysis of the

341

transcriptomic profilings of Mandarin fish (Siniperca chuatsi) infected with

342

Flavobacterium columnare with an emphasis on immune responses. Fish

343

Shellfish immun. 2015; 43(1):111-9.

M AN U

SC

340

344

[11] Wang Y, Lu Y, Zhang Y, Ning Z, Li Y, Zhao Q, et al. The draft genome of the

345

grass carp (Ctenopharyngodon idellus) provides insights into its evolution and

346

vegetarian adaptation. Nat genet. 2015; 47(6):625-31.

349 350

TE D

348

[12] Kim D, Langmead B, Salzberg S L. HISAT: a fast spliced aligner with low memory requirements. Nat methods. 2015; 12(4):357-60. [13] Anders S, Pyl P T, Huber W. HTSeq—a Python framework to work with high-throughput sequencing data. Bioinformatics. 2014:btu638.

EP

347

[14] Wang L, Feng Z, Wang X, Wang X, Zhang X. DEGseq: an R package for

352

identifying differentially expressed genes from RNA-seq data. Bioinformatics,

353

2010; 26(1):136-38.

354 355

AC C

351

[15] Jia Y, Lin F, Hongkun Z, Yong Z, Jie C, Zengjin Z, et al. WEGO: a web tool for plotting GO annotations. Nucleic Acids Res. 2006; 34(2):: 293-7.

356

[16] Harris MA, Clark J, Ireland A, Lomax J, Ashburner M, Foulger R, et al. The

357

Gene Ontology (GO) database and informatics resource. Nucleic Acids Res.

358

2004; 32(suppl 1):258-61.

359

[17] Minoru K, Susumu G, Miho F, Mao T, Mika H. KEGG for representation and

360

analysis of molecular networks involving diseases and drugs. Nucleic Acids Res.

361

2010; 38(1):355-60.

ACCEPTED MANUSCRIPT 362

[18] Litman GW, Yoder JA, Cannon JP, Haire RN. Novel Immune-type Receptor

363

Genes and the Origins of Adaptive and Innate Immune Recognition1. Integr

364

Comp Biol. 2009; 41(6):1507-8. [19] Sato S, Miller AS, Inaoki M, Bock CB, Jansen PJ, Tang MLK, et al. CD22 Is

366

Both a Positive and Negative Regulator of B Lymphocyte Antigen Receptor

367

Signal Transduction: Altered Signaling in CD22-Deficient Mice. Immunity.

368

1996; 5(6):551-62.

371 372 373 374 375 376

response. Central-Euro J Immun. 2014; 1(1):109-15.

SC

370

[20] Osińska I, Popko K, Demkow U. Perforin: an important player in immune

[21] Aderem A, Underhill DM. Mechanisms of phagocytosis in macrophages. Annu Rev Immunol. 1999; 17593-623.

M AN U

369

RI PT

365

[22] Jutras I, Desjardins M. Phagocytosis: at the crossroads of innate and adaptive immunity. Annu Rev Cell Dev Bio. 2005; 21(21):511-27. [23] Martinez FO, Sica AA, Locati M. Macrophage activation and polarization. Front Biosci. 2008; 13(4):453-61.

[24] Vandivier RW, Ogden CA, Fadok VA, Hoffmann PR, Brown KK, Botto M, et al.

378

Role of surfactant proteins A, D, and C1q in the clearance of apoptotic cells in

379

vivo and in vitro: calreticulin and CD91 as a common collectin receptor

380

complex. J Immunol. 2002; 169(7):3978-86.

383 384

EP

382

[25] Hynes RO. Integrins: versatility, modulation, and signaling in cell adhesion. Cell. 1992; 69(1):11-25.

[26] Hampton MB, Kettle AJ, Winterbourn CC. Inside the neutrophil phagosome:

AC C

381

TE D

377

oxidants, myeloperoxidase, and bacterial killing. Blood. 1998; 92(9):3007-17.

385

[27] Zhi JX, Yuan L, yu MW. Research advance of the lysozyme. Lett Biotech. 2002.

386

[28] Repnik U, Stoka V, Turk V, Turk B. Lysosomes and lysosomal cathepsins in cell

387

death. BBA Proteins Proteomics. 2012; 1824(1):22-33.

388

[29] Turk V, Stoka V, Vasiljeva O, Renko M, Sun T, Turk B, et al. Cysteine

389

cathepsins: From structure, function and regulation to new frontiers. Biochimica

390

Et Biophysica Acta. 2012; 1824(1):68-88.

391

[30] He J, Tohyama Y, Yamamoto KI, Kobayashi M, Shi Y, Takano T, et al.

ACCEPTED MANUSCRIPT 392

Lysosome is a primary organelle in B cell receptor-mediated apoptosis: an

393

indispensable role of Syk in lysosomal function. Genes Cells. 2005;

394

10(1):23-35.

397 398 399 400

immun. 2002; 12(5):399-420.

RI PT

396

[31] Holland MC, Lambris JD. The complement system in teleosts. Fish Shellfish

[32] Boshra H, Li J, Sunyer J. Recent advances on the complement system of teleost fish. Fish Shellfish immun. 2006; 20(2):239-62.

[33] Ricklin D, Lambris JD. Complement in immune and inflammatory disorders: pathophysiological mechanisms. Journal of Immunology. 2013; 190(8):3831-8.

SC

395

[34] Gando S, Nanzaki S, Sasaki S, Aoi K, Kemmotsu O. Activation of the extrinsic

402

coagulation pathway in patients with severe sepsis and septic shock. Crit Care

403

Med. 1998; 26(12):2005-9.

404 405

M AN U

401

[35] Lizbeth H. Serine protease mechanism and specificity. Chem Rev. 2002; 102(6):4501-23.

[36] Orfeo T, Brummel-Ziedins KE, Butenas S, Mann KG. Simulating blood

407

coagulation: the contribution of α2-macroglobulin and α1-antitrypsin. Faseb J.

408

2006; 20(4):A64.

TE D

406

[37] Song X, Zhao J, Bo Y, Liu Z, Wu K, Gong C. Aeromonas hydrophila induces

410

intestinal inflammation in grass carp ( Ctenopharyngodon idella ): An

411

experimental model. Aquaculture. 2014; 434171-8.

413 414 415 416 417 418 419 420 421

[38] Salazar-Mather TP, Hokeness KL. Cytokine and Chemokine Networks: Pathways to Antiviral Defense. Curr Top Microbio. 2006; 303(1):29-46.

AC C

412

EP

409

[39] Savan R, Sakai M. Genomics of fish cytokines. Comp Biochem Phys D. 2006; 1(1):89-101.

[40] Robertsen B. The interferon system of teleost fish. Fish Shellfish immun. 2006; 20(2):172-91. [41] Secombes CJ, Wang T, Bird S. The interleukins of fish. Dev Comp Immunol. 2011; 35(12):1336-45. [42] Alejo A, Tafalla C. Chemokines in teleost fish species. Dev Comp Immunol. 2011; 35(12):1215-22.

ACCEPTED MANUSCRIPT 422

[43] Su J, Heng J, Teng H, Peng L, Yang C, Li Q. Identification, mRNA expression

423

and genomic structure of TLR22 and its association with GCRV

424

susceptibility/resistance in grass carp (Ctenopharyngodon idella). Dev Comp

425

Immunol. 2011; 36(2):450-62. [44] Hu GB, Zhang SF, Xi Y, Liu DH, Liu QM, Zhang SC. Cloning and expression

427

analysis of a Toll-like receptor 22 ( tlr22 ) gene from turbot, Scophthalmus

428

maximus. Fish Shellfish immun. 2015; 44(2):399-409.

429

SC

430 431

M AN U

432 433 434 435 436

441 442 443 444 445 446 447 448 449 450

EP

440

AC C

439

TE D

437 438

RI PT

426

ACCEPTED MANUSCRIPT Figure legends

452

Fig. 1. (A. B) Overview of differently expressed genes distribution trends between

453

infected (YC_0909) and non-infected (YD_0909) group. (A) Gene expression

454

distribution. The x-axis is gene expression quantities (in logs) of non-infected group

455

and y-axis is gene expression quantities of infected group. (B) Volcano Plot. The log2

456

fold (infected/ non-infected) is indicated the mean expression level for each genes.

457

Each dot represents one gene. Red dots represent up-regulated genes and blue dots

458

represent down-regulated genes. Green dots represent no differential expression

459

genes.

460

Fig. 2. Gene ontology (GO) enrichment analysis of differently expressed genes. The

461

results of GO enrichment analysis of differently expressed genes were classified into

462

three categories: biological process, cellular component and molecular function. The

463

x-axis is gene functional classification of GO, y-axis is the corresponding number of

464

up-regulated and down-regulated genes.

465

Fig. 3. The top30 most enriched KEEG pathways. The x-axis is KEGG pathway

466

classification and y-axis is the corresponding number of up-regulated and

467

down-regulated genes.

468

Fig. 4. Comparison of the expressions of RNA-Seq and qRT-PCR results. The

469

transcript expression levels of the selected genes were each normalized to that of the

470

β-actin gene.

AC C

EP

TE D

M AN U

SC

RI PT

451

ACCEPTED MANUSCRIPT Table 1 Primers used for qRT-PCR verification of differently expressed genes

TGATGGCTTTCTCAC ACTGC CTTGTACCGAGTCGG IL-1β ATGGT GCAACTGGGCTCAA TNF-α GCTTAC TCCTCGTGTGAAAGC LysC AAGAC GTGGAGCCGAATGT F7 CCTAAA GTTTCATGCATCGTC C1 CAAAG GCCATGGCCAGTAA C3 CTTTGT GGCAAAGGGTCAGT C7 GTGTTT GGCTTCTATCAGCTG CD 22 CTGTG TCCTTACGAGGCATG TLR 22 AGCTT ATGACTTCGATGAG IRF 4 CTGGTG TTCAAGCTCCGTTAA MPEG 1 TGCTG ATCACCTGGCAGAA MHC-I Ag AAATGG TCTTCCCTCCTCCTG MHC-II Ag TCAGA AAATCCAGAAGCTG HSP 70 CTCCAA GGCCTTCTGTCTGTC α2MG CTCTG GTGCCCATCTACGA β-actin GGGTTA

GATATTTCTTCGCCGTT GGA GAACAGGAGGTTGGCA TGTT GGTCCTGGTTCACTCT CCAA ATCCCTCAAATCCATC AAGC GTCATGCTCACCTGCG ACTA TAAACAGCCAGTCCAC CAGA CTCCATGAACGACAGC TGAA ATTTGGTTTCTGGCCA ACTG CCACCGAAGATCCTTC AAAT CGACAAGAGGAGGGT GAGAG CGTGTGTCGTGTTTTCT TCC TTTCACTTGGTGGGTCT TCA CCTGATGTTCCACCAC ACAG CCCTCCACAGGTGTGA ACTT TCAATTCCAAGGGACA GAGG AACCATGATGCACTTG GACA TCTCAGCTGTGGTGGT GAAG

AC C

EP

TE D

M AN U

XCR 1

Reverse primer (5′-3′)

Accession number KF937391.1 JQ692172.1

RI PT

Forward primer (5′-3′)

SC

Gene

JQ670915.1

AF402599.1

JX088706.1 JQ358795.1 AY374472.1 JN655169.1 NM_0010838 23.2 HQ676542.2 NM_0011227 10.1 NM_212737. 1 AB221535.1 EF140725.1 EU816595.1 AY425709.1 M25013.1

ACCEPTED MANUSCRIPT Table 2 Summary of sequencing and mapping Non-infected Group

Read Length (bp)

125

125

Total Raw Reads

65513172

63438956

Total Clean Reads

64627432

62568036

Mapped reads Ratio (%)

79.88

81.41

Detected Gene Number

24507

24383

AC C

EP

TE D

M AN U

SC

RI PT

Infected Group

ACCEPTED MANUSCRIPT Table 3 List of the differently expressed genes (DEGs) involved in immune responses Category and Gene name

Log2Fold

Diff

Qvalue

1.807160033

Up

5.12E-276

Calmegin precursor

1.72725785

Up

8.66E-35

Vesicle-trafficking protein SEC22b-A

1.66135348

Up

8.89E-38

C-type lectin domain family 4 member

1.578698141

Up

8.97E-143

Calreticulin precursor

1.452254815

Up

0

C-type lysozyme

1.424010322

Up

1.65E-20

Up

0

Up

4.22E-12

Up

9.89E-23

Transferrin receptor 1a

Tubulin beta-1 chain

1.38832028

RI PT

Phagocytosis

1.344612409

Integrin beta 3a

1.330257571

Tubulin alpha-1C chain

1.283833515

Up

3.85E-08

Macrophage expressed gene 1 (MPEG 1)

1.237597456

Up

1.90E-43

Cysteine proteinase

1.217510291

Up

2.12E-40

1.198804756

Up

1.80E-97

1.13560686

Up

4.78E-40

M AN U

SC

Macrophage mannose receptor 1

Cathepsin K Neutrophil cytosolic factor 1 Cathepsin L precursor

1.12855511

Up

0.000121212

1.088007207

Up

1.15E-68

1.074417256

Up

5.08E-233

1.073145635

Up

4.97E-180

1.06144093

Up

7.49E-08

1.026675458

Up

7.77E-07

Integrin alpha-V precursor

1.011633801

Up

2.12E-65

Neutrophil cytosolic factor 4

1.008973284

Up

0.000100343

Complement Component C3b

4.309367861

Up

8.24E-05

Factor B (Bf)

2.298309327

Up

0

Complement component C1

2.168355572

Up

5.50E-07

Complement component C7

1.438283815

Up

0

1.346203375

Up

0

Coagulation factor Ⅶ

-1.07107892

Down

1.77E-210

Coagulation factor Ⅸ precursor

-1.07927431

Down

2.53E-48

G-type lysozyme Tubulin alpha-1B chain Cathepsin B, precursor Coronin-2B isoform X2

EP

Complement system

TE D

Cystinosin precursor

AC C

Complement component C8β

precursor

Complement component 3 (C3)

-1.11192977

Down

4.05E-63

Mannose-binding serine protease

-1.17167819

Down

6.27E-28

Small inducible cytokine B14 precursor

-1.43605885

Down

0.000921317

α-2-macroglobulin

-1.49752246

Down

0

Complement factor H like 4 precursor

-3.11609747

Down

7.40E-60

Tumor necrosis factor-α induced protein 6

5.346842986

Up

5.49E-18

C-X-C motif chemokine 13 precursor

2.280377621

Up

3.44E-227

Interferon regulatory factor 4 (IRF4)

2.06144081

Up

1.91E-09

1.641029822

Up

1.41E-121

Cytokine

Chemokine-like receptor 1

ACCEPTED MANUSCRIPT Category and Gene name

Log2Fold

Diff

Qvalue

1.596217893

Up

5.32E-09

Perforin-like protein 2

1.415077953

Up

4.55E-05

Interferon-induced transmembrane protein 1

1.349685778

Up

4.56E-15

Matrix metalloproteinase 9

1.279235669

Up

0

Chemokine C motif receptor XCR1

1.204045397

Up

7.61E-09

Interleukin-21 receptor precursor

1.123368854

Up

6.71E-06

Interleukin-10

1.114121126

Up

2.99E-11

C-X-C motif chemokine 9-like

1.113733102

Up

3.67E-191

Granulocyte colony-stimulating factor receptor

1.064605775

Up

1.52E-53

Tumor necrosis factor receptor 2

1.058231015

Up

6.22E-18

Interferon regulatory factor 10

1.057900539

Up

2.48E-16

Interleukin-1 receptor I

1.055859594

Up

2.10E-20

CXCR3

1.001675808

Up

7.29E-24

Interleukin-27 subunit beta precursor

-1.35965669

Down

2.17E-100

Interleukin 13 receptor, alpha 1 precursor

-1.52185706

Down

1.42E-56

-2.02602217

Down

0.000352238

4.148903699

Up

0.00024541

3.391288769

Up

1.26E-192

3.061440646

Up

1.38E-13

2.338280491

Up

2.13E-08

1.769272742

Up

0

1.687275774

Up

1.06E-28

1.491680658

Up

3.26E-11

Novel immune-type receptor 5 (NITR5)

1.470246033

Up

8.15E-06

Perforin-like protein 2

1.415077953

Up

4.54746E-05

1.397043842

Up

3.00E-10

1.361001026

Up

2.42E-06

1.313827508

Up

3.06E-18

1.243644382

Up

0.000768153

CD180 antigen-like

1.208281846

Up

0.000747575

CD2 family receptor (CD2 f)

1.119543937

Up

3.11E-07

CD209 antigen

1.077742201

Up

7.19E-08

Lymphocyte activation gene 3 precursor

1.049805976

Up

2.13E-19

1.04181149

Up

1.07E-11

proteasome activator complex subunit 2

1.006913826

Up

2.92E-48

CD36 antigen

-1.01129673

Down

1.20E-90

MHC class II alpha antigen

-1.51012641

Down

2.54E-103

MHC class II DA-beta-2, partial

-1.60152434

Down

7.05E-51

3.203892677

Up

Log2Fold

Diff

SC

M AN U

Immunity-related GTPase family, LRG-47 Antigen processing and presentation MHC class I antigen Lymphocyte antigen 6F-like T cell receptor beta chain Novel immune-type receptor protein CD79a precursor B-cell receptor CD22

CD11

EP

Perforin 3

TE D

Heat shock protein 70 (HSP70)

Lymphocyte antigen 86 precursor

T-lymphoma invasion and metastasis-inducing

AC C

protein

CD9 antigen-like

RI PT

Interleukin-1 beta

Pattern recognition receptor Intelectin Category and Gene name

0 Qvalue

ACCEPTED MANUSCRIPT C-type lectin 4

2.39295431

Up

2.41E-19

Toll-like receptor 22

1.70110923

Up

9.10E-18

NALPL2

1.154550248

Up

1.10E-10

C-type lectin 10

1.118247028

Up

2.11E-22

NLRC3 protein

1.099242548

Up

2.90E-10

JAM-C-like protein, partial

3.231365769

Up

0.000613005

T cell receptor beta chain, partial

3.061440646

Up

1.38E-13

Claudin-1 precursor

2.563940973

Up

2.01E-05

Versican b precursor

1.831958763

Up

3.29E-31

Up

5.48E-05

Up

1.30E-13

Leukocyte immune-type receptor

1.38336861

RI PT

Cell adhesion molecules

Selectin P ligand precursor

1.119156551

Claudin b, partial

-1.46580617

Claudin 2

-1.65737772

Down

8.21E-24

Novel protein similar to E-cadherin (cdh1)

-2.23073992

Down

3.01E-235

Claudin-4

-3.63899852

Down

1.69E-05

1.604233798

Up

5.05E-35

1.440956114

Up

2.14E-36

1.365951445

Up

6.22E-15

1.240224556

Up

4.98E-51

Oxidoreductase-like domain-containing protein

1.616029493

Up

0.000228787

Peroxisomal membrane protein 4

1.588688079

Up

6.95E-08

Glutathione synthetase

1.567793403

Up

7.30E-36

1.146102839

Up

3.10E-09

1.13170309

Up

4.08E-49

1.012681682

Up

1.67E-23

-1.0320958

Down

3.73E-56

-1.30329752

Down

0

-1.37843725

Down

2.83E-83

Reticulocalbin-3 precursor Caspase recruitment domain protein Caspase-7

EF-hand calcium-binding domain-containing protein

Spermidine synthase

TE D

Antioxidant enzymes

Natural killer cell enhancing factor Arginase 2C precursor Catalase

EP

Glutathione S-transferase

Gamma-glutamylcyclotransferase precursor

AC C

SC

M AN U

Apoptosis

Down

0.000656937

AC C

EP

TE D

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

AC C

EP

TE D

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

AC C

EP

TE D

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

AC C

EP

TE D

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

ACCEPTED MANUSCRIPT Highlights 1. Transcriptome analysis of the grass carp infected A. hydrophila was performed. 2. 2121 differential expression genes were identified. 3. 105 important immune related genes were categorized and discussed.

AC C

EP

TE D

M AN U

SC

RI PT

4. For a better understanding of immune mechanism of grass carp infected bacteria.