Abstract

   The Japanese scallop (Mizuhopecten yessoensis) is one of the most
   important aquaculture species in Asian countries; however, it has
   suffered severe infection by Polydora in northern China in recent
   years, causing great economic losses. The Polydora parasitizes the
   shell of scallops, badly destroying the shell’s structure. To
   investigate the molecular response mechanism of M. yessoensis to
   Polydora infestion, a comprehensive and niche-targeted cDNA sequence
   database for diseased scallops was constructed. Additionally, the
   transcriptional changes in the edge mantle, central mantle and
   hemocytes, tissues directly related to the disease, were first
   described in this study. The results showed that genes involved in
   shell formation and immunomodulation were significantly differentially
   expressed due to the infestation. Different transcriptional changes
   existed between the two mantle regions, indicating the different
   molecular functions likely responsible for the formation of different
   shell layers. The differential expression of genes for immune
   recognition, signal transduction and pathogen elimination presented an
   integrated immune response process in scallops. Moreover,
   neuromodulation and glycometabolism involved in the regulation process
   with relevant function significantly enriched. The study provides
   valuable information for mechanism study of shell formation and
   immunomodulation in scallops.

Introduction

   The Japanese scallop (Mizuhopecten yessoensis) is a large and old
   (dating back to ~350 Ma) group living on the cold and stable ocean
   bottoms of the northwestern Pacific Ocean. This species has been a
   major economic aquaculture species in Asian countries and is consumed
   worldwide^[36]1. However, the frequent outbreak of bacterial and
   parasitic diseases in M. yessoensis has caused severe economic losses.
   Spionid worms of the genus Polydora, which mainly inhabit the shells of
   many benthic bivalves^[37]2–[38]5 including some commercially important
   species, have become a serious problem in recent scallop aquaculture.
   It was roughly estimated that over 90% of cultured M. yessoensis in the
   Dalian Zhangzidao Sea (Liaoning, China) had been infested by Polydora
   at various degrees, badly affecting the development of the scallop
   industry. Polydora mainly parasitize the left valve of the scallop by
   excavating U-type tubes, which badly damages the structure of the
   shells. To prevent further invasion, scallops accelerate the secretion
   of shell components to repair the damaged regions. As a result, many
   protuberances are formed on the inner surface of the shell
   (Fig. [39]1), which are signs of the shell repair process. When
   Polydora drill through the shell, they can directly infect the soft
   body (e.g. the mantle tissue) and the damaged shell exposes the scallop
   to various pathogens, activating the scallop’s immune modulation to
   combat the infection. At present, research about the Polydora disease
   mainly focuses on the biological characteristics and parasitic
   behaviour of Polydora^[40]6–[41]11; however, the response of M.
   yessoensis, especially on the molecular level, has never been reported,
   limiting the understanding of disease resistance in scallops.

Figure 1.

   Figure 1
   [42]Open in a new tab

   The outer and inner shell surfaces of the healthy (A) and infested by
   Polydora (B) M. yessoensis.

   The mantle is the main tissue involved in the shell formation, located
   midway between the shell and visceral mass. Different mantle regions
   are responsible for the formation of different shell layers with the
   edge mantle being responsible for the periostracal and prismatic
   layerand the central mantle for the nacreous layer^[43]12–[44]14. The
   resistance of M. yessoensis to the Polydora disease is characterized by
   the shell repair process with protuberance formation on the inner
   surface. Mantle tissue is therefore essential in the defense against
   Polydora. In recent years, the molecular function of mantle tissue has
   been extensively studied. With the prevalence of next-generation
   sequencing, transcriptomes of the mantle tissue have been sequenced in
   many mollusks, such as Patella vulgate^[45]15, Mytilus edulis^[46]16,
   Hyriopsis cumingii^[47]17, Pinctada fucata^[48]18, P.
   margaritifera^[49]19,[50]20, Chlamys farreri^[51]21, and M.
   yessoensis^[52]22,[53]23, to explore the mechanisms of shell formation.
   A wide variety of candidate genes involved in biomineralization have
   been identified, such as calcium-binding protein, carbohydrate-binding
   protein, glycoprotein, collagen, nacrein, perlucin and others, and
   their functions in shell formation have been partially
   characterized^[54]15–[55]23. However, due to the extreme complexity of
   shell formation and the differences among species, our understanding is
   still fragmented and lacks consensus. In addition, transcriptional
   changes in the mantle tissue of Ruditapes philippinarum caused by Brown
   Ring Disease (a shell disease caused by the bacteria Vibrio tapetis)
   not only include genes related to biomineralization but also many
   immunity genes^[56]24,[57]25 (such as big defensin, lectin, C1q domain
   containing protein and so on), implying the important role of mantle
   tissue in the immune response to pathogens.

   Like all invertebrates, scallops rely on an exclusively innate immune
   system to execute cellular and humoural (macromolecules) immune
   reactions to resist invaders^[58]26. As the classical immune tissue in
   scallops, hemocytes use encapsulation and phagocytosis to eliminate
   undesirable particles (e.g., pathogens and abnormal/dead cells)
   performing the cellular immune reactions and contribute to the
   secretion of immune macromolecules to execute humoural immune
   responses, such as immune recognition, signal transduction and
   destruction of invaders, as reviewed by Song et al.^[59]27. Recent
   studies have also demonstrated a key role of hemocytes in
   biomineralization due to their involvement in mineral transport, as
   well as production of the extracellular matrix of the
   shells^[60]28–[61]31, indicating the multifunction of hemocytes.
   Ivanina et al., recognized that potential trade-offs existed between
   biomineralization and immunity in Crassostrea gigas and C.
   virginica^[62]29, but relevant studies have not been reported in M.
   yessoensis.

   In the present study, a comprehensive and niche-targeted cDNA sequences
   database was constructed for M. yessoensis infested by Polydora,
   transcriptional changes in the edge mantle, central mantle and
   hemocytes of the diseased scallops were first detected, and genes
   related to shell formation and immunity were screened. The aim was to
   reveal the molecular response mechanism of M. yessoensis to the
   infestation of Polydora, and provide more valuable information for the
   mechanism studies of biomineralization and immunomodulation.

Results

Transcriptome sequencing, assembly and functional annotation

   To identify differently expressed genes between healthy and diseased M.
   yessoensis infested by Polydora, 18 cDNA libraries of the edge mantle,
   central mantle and hemocytes of two groups were separately constructed
   for RNA-seq sequencing. After removing adaptors and low-quality reads,
   a total of 942,738,120 reads (~109.7 Gb) with an average of 52,374,340
   (~6.1 Gb) for each sample were obtained, and the detailed sequencing
   information for each library is listed in Supplementary Table [63]S1.
   De novo assembly for the transcriptome sequences was chosen to
   construct a comprehensive and niche-targeted cDNA sequence database for
   diseased scallops, and the Trinity method was used with all trimmed
   reads. Finally, 80,831 unigenes were obtained with an average length of
   1607 bp (Supplementary Fig. [64]S1), and the assembly information is
   summarized in Supplementary Table [65]S2. The unigenes were annotated
   by searching the sequences against the Nr, SWISS-PROT, KOG, GO and KEGG
   databases using BLASTX with a cut-off of Evalue ≤ 1e-5 to describe
   their functions at different levels, and the annotation ratios were
   31.32%, 23.75%, 20.27%, 21.73% and 8.97%, respectively (Supplementary
   Table [66]S3).

   For GO classification, 17,567 unigenes were annotated into 12,658 GO
   terms in three main GO categories: biological process (~65.67%),
   cellular component (~10.44%), and molecular function (~23.90%)
   (Supplementary Table [67]S4 and Fig. [68]S2). In the biological process
   category, most of the unigenes were assigned to the terms of ‘cellular
   process’ (~14.81%), ‘single-organism process’ (~12.65%) and ‘metabolic
   process’ (~11.58%), including ‘immune system process’ (~1.15%). In the
   cellular component category, most of the unigenes were assigned to the
   terms of ‘cell’ (~19.54%), ‘cell part’ (~19.49%) and ‘organelle part’
   (~15.29%). In the molecular function category, most of the unigenes
   were assigned to the terms of ‘binding’ (~44.48%), ‘catalytic activity’
   (~33.17%) and ‘transporter activity’ (~5.95%). A total of 16,387
   unigenes were assigned to 25 ortholog groups in the KOG database
   (Supplementary Table [69]S5 and Fig. [70]S3), which were distributed
   mainly into three groups: (R) ‘General function prediction only’
   (~19.30%), (T) ‘Signal transduction mechanisms’ (~17.82%), and (O)
   ‘Posttranslational modification, protein turnover, chaperones’
   (~10.02%). The biological signaling pathways for the unigenes were
   predicted by the KEGG database, and 7253 unigenes were assigned to 342
   pathways (Table [71]S6) involving six KEGG categories with ‘cellular
   processes’ (~11.81%), ‘environmental information processing’ (~13.76%),
   ‘genetic information processing’ (~7.56%), ‘human diseases’ (~25.05%),
   ‘metabolism’(~18.38%) and ‘organismal systems’ (~23.43%), in which the
   subcategories of ‘infectious diseases’, ‘immune diseases’, ‘nervous
   system’, ‘immune system’, etc. were contained (Supplementary
   Fig. [72]S4).

Overall gene regulation in diseased scallops

   The expression levels of the assembled unigenes were estimated by the
   FPKM method. Gene differential expression analysis was conducted
   between healthy and infested scallops in the edge mantle, central
   mantle and hemocytes tissues (p ≤ 0.05). A total of 3670 unigenes were
   differentially expressed in at least one of the three tissues with 1729
   in the edge mantle, 1992 in the central mantle and 1010 in the
   hemocytes. Among these DEGs, there were 822 DEGs (~22.4%) shared by at
   least two of the tissues with 623 between the edge and central mantles,
   314 between the edge mantle and hemocytes and 363 between the central
   mantle and hemocytes. Moreover, 239 (~6.5%) DEGs were shared among all
   the three tissues (Fig. [73]2). It is noteworthy that these
   three-tissue-shared DEGs were significantly enriched in the functions
   and processes related to glycometabolism by GO and KEGG enrichment
   analysis (p ≤ 0.05) (Figs [74]3 and [75]4). For instance, the GO terms
   of the ‘fructose metabolic process’, ‘glucose 6−phosphate metabolic
   process’, ‘cellular glucose homeostasis’, ‘mannose metabolic process’,
   ‘glucose metabolic process’ and ‘glycolytic process’ were the most
   significantly enriched terms in the biological process category
   (Fig. [76]3A). Meanwhile, relevant molecular activity, such as
   ‘mannokinase activity’, ‘hexokinase activity’, ‘fructokinase activity’,
   ‘glucokinase activity’, ‘glucose binding’ and ‘ATP binding’ were also
   significantly overrepresented in the molecular function category
   (Fig. [77]3B). Additionally, KEGG enrichment analysis showed a similar
   trend with pathways related to glycometabolism being the most
   significantly enriched, including ‘carbohydrate digestion and
   absorption’, ‘galactose metabolism’, ‘fructose and mannose metabolism’,
   ‘starch and sucrose metabolism’, ‘glycolysis/gluconeogenesis’ and
   ‘carbon metabolism’ (Fig. [78]4). Moreover, all of the above relevant
   terms were found to be upregulated in all the three tissues
   (Tables [79]1–[80]3), which likely implied very active energy
   metabolism in these tissues to meet the living demands of the diseased
   scallops.

Figure 2.

   Figure 2
   [81]Open in a new tab

   Distribution of differentially expressed unigenes between diseased and
   healthy M. yessoensis across different tissues (edge mantle, central
   mantle and hemocytes).

Figure 3.

   [82]Figure 3
   [83]Open in a new tab

   GO terms significantly enriched among DEGs shared by the edge mantle,
   central mantle and hemocytes. (A) Biological process. (B) Molecular
   function. Only the top 20 are presented, which are sorted by the
   p-value (p ≤ 0.05) and displayed with the –log[10](p-value). Terms or
   pathways related to ‘glycometabolism’ and ‘immune response’ are
   separately highlighted with green and orange dots.

Figure 4.

   Figure 4
   [84]Open in a new tab

   KEGG pathway enrichment analysis among DEGs shared by the three
   tissues. Only the top 20 are presented, which are sorted by the p-value
   (p ≤ 0.05). Terms or pathways related to ‘glycometabolism metabolism’
   are highlighted with green triangles.

Table 1.

   GO terms significantly enriched in the edge mantle.
   iD Term p(≤0.05) Type iD Term p(≤0.05) Type
   Upregulated in diseased scallops Related to cellular motility
   Related to calcium ion GO:0005929 cilium 2.29E-12 C
   GO:0060402 calcium ion transport into cytosol 6.53E-04 P GO:0005874
   microtubule 8.31E-11 C
   GO:0005509 calcium ion binding 8.02E-04 F GO:0031514 motile cilium
   1.49E-10 C
   GO:0051480 cytosolic calcium ion homeostasis 2.81E-03 P GO:0003777
   microtubule motor activity 2.89E-08 F
   GO:0005262 calcium channel activity 5.00E-02 F GO:0005930 axoneme
   4.39E-08 C
   Related to shell organic matrix GO:0003341 cilium movement 8.41E-08 P
   GO:0018392 glycoprotein 3-alpha-L-fucosyltransferase activity 1.88E-03
   F GO:1990716 axonemal central apparatus 1.08E-07 C
   GO:0030203 glycosaminoglycan metabolic process 1.40E-02 P GO:0035082
   axoneme assembly 6.89E-07 P
   Related to nervous system GO:0030286 dynein complex 1.42E-06 C
   GO:0021591 ventricular system development 4.77E-05 P GO:0060294 cilium
   movement involved in cell motility 1.27E-05 P
   GO:0021670 lateral ventricle development 1.15E-04 P GO:0003351
   epithelial cilium movement 3.22E-05 P
   GO:0043083 synaptic cleft 4.64E-04 C GO:0036158 outer dynein arm
   assembly 3.22E-05 P
   GO:0007399 nervous system development 4.75E-04 P GO:0036157 outer
   dynein arm 5.83E-05 C
   GO:0010976 positive regulation of neuron projection development
   2.96E-03 P GO:0010632 regulation of epithelial cell migration 6.79E-05
   P
   GO:0007420 brain development 4.64E-03 P GO:0007018 microtubule-based
   movement 9.43E-05 P
   GO:0007417 central nervous system development 4.68E-03 P GO:0060271
   cilium morphogenesis 4.77E-04 P
   GO:0006665 sphingolipid metabolic process 7.02E-03 P GO:0036064 ciliary
   basal body 6.35E-04 C
   GO:0045202 synapse 1.76E-02 C GO:0042384 cilium assembly 1.09E-03 P
   GO:0043005 neuron projection 3.21E-02 C GO:0036156 inner dynein arm
   1.50E-03 C
   GO:0005328 neurotransmitter:sodium symporter activity 1.39E-02 F
   GO:0001539 cilium or flagellum-dependent cell motility 2.32E-03 P
   Related to energy metabolism GO:0005858 axonemal dynein complex
   2.81E-03 C
   GO:0004550 nucleoside diphosphate kinase activity 3.94E-08 F GO:0036159
   inner dynein arm assembly 2.81E-03 P
   GO:0006165 nucleoside diphosphate phosphorylation 3.62E-06 P GO:0008017
   microtubule binding 1.54E-02 F
   GO:0009142 nucleoside triphosphate biosynthetic process 2.08E-05 P
   GO:0072372 primary cilium 3.60E-02 C
   GO:0016887 ATPase activity 2.58E-03 F GO:0030334 regulation of cell
   migration 3.84E-02 P
   GO:0004396 hexokinase activity 1.15E-04 F
   GO:0019158 mannokinase activity 1.15E-04 F Downregulated in diseased
   scallops
   GO:0004340 glucokinase activity 3.13E-04 F Related to nervous system
   GO:0006000 fructose metabolic process 3.13E-04 P GO:0043523 regulation
   of neuron apoptotic process 2.40E-03 P
   GO:0008865 fructokinase activity 3.13E-04 F Related to immume system
   GO:0016051 carbohydrate biosynthetic process 5.06E-04 P GO:0006954
   inflammatory response 2.07E-03 P
   GO:0005536 glucose binding 6.53E-04 F GO:0051607 defense response to
   virus 6.20E-03 P
   GO:0001678 cellular glucose homeostasis 8.86E-04 P GO:0045087 innate
   immune response 1.00E-02 P
   GO:0051156 glucose 6-phosphate metabolic process 8.86E-04 P GO:0005044
   scavenger receptor activity 1.31E-02 F
   GO:0006013 mannose metabolic process 4.67E-03 P
   GO:0006006 glucose metabolic process 2.64E-02 P
   [85]Open in a new tab

Table 3.

   GO terms significantly enriched in the hemocytes.
   iD Term p(≤0.05) Type iD Term p(≤0.05) Type
   Upregulated in diseased scallops Downregulated in diseased scallops
   Related to energy metabolism Related to nervous system
   GO:0004396 hexokinase activity 4.46E-06 F GO:0007218 neuropeptide
   signaling pathway 5.24E-03 P
   GO:0019158 mannokinase activity 4.46E-06 F GO:0008188 neuropeptide
   receptor activity 9.90E-03 F
   GO:0004340 glucokinase activity 1.24E-05 F Related to immume response
   GO:0006000 fructose metabolic process 1.24E-05 P GO:0002218 activation
   of innate immune response 3.53E-05 P
   GO:0008865 fructokinase activity 1.24E-05 F GO:0051607 defense response
   to virus 1.24E-03 P
   GO:0005536 glucose binding 2.63E-05 F GO:0002230 positive regulation of
   defense response to virus by host 1.75E-03 P
   GO:0001678 cellular glucose homeostasis 3.60E-05 P GO:0045087 innate
   immune response 1.52E-02 P
   GO:0051156 glucose 6-phosphate metabolic process 3.60E-05 P GO:0005044
   scavenger receptor activity 1.67E-02 F
   GO:0006013 mannose metabolic process 2.04E-04 P GO:0006954 inflammatory
   response 1.70E-02 P
   GO:0006006 glucose metabolic process 6.57E-03 P GO:0006955 immune
   response 3.66E-02 P
   GO:0006096 glycolytic process 1.19E-02 P
   GO:0016887 ATPase activity 1.67E-02 F
   GO:0005524 ATP binding 2.11E-02 F
   Related to nervous system
   GO:0005328 neurotransmitter:sodium symporter activity 1.23E-02 F
   [86]Open in a new tab

Gene regulation in different tissues

   For the DEGs of the three tissues, there were 1162 (~67.21%)
   upregulated and 567 (~32.79%) downregulated in the edge mantle, 1040
   (~52.21%) upregulated and 952 (~47.79%) downregulated in the central
   mantle, and 521 (~51.58%) upregulated and 489 (~48.42%) downregulated
   in hemocytes. GO enrichment analysis was separately carried out for the
   up and downregulated DEGs of each tissue. On the whole, a variety of
   functions were covered for the enriched GO terms, but many of them
   could be classified into one category according to their functions
   (Tables [87]1–[88]3, Supplementary Fig. [89]5S), which were probably
   related to shell formation and immunomodulation, such as calcium ion
   binding and transport, shell organic matrix, nervous system
   development, epithelial cell motility and immune response, including
   the glycometabolism-related functions shared by the three tissues.

   GO enrichment analysis for upregulated genes in the edge mantle are
   described in Table [90]1. First, some GO terms related to calcium ion
   activities, such as ‘calcium ion transport into cytosol’, ‘calcium ion
   binding’, ‘cytosolic calcium ion homeostasis’ and ‘calcium channel
   activity’ and shell organic matrix, such as ‘glycoprotein
   3-alpha-L-fucosyltransferase activity’ and ‘glycosaminoglycan metabolic
   process’, were significantly enriched, which were closely related to
   the process of biomineralization. In addition, many GO terms related to
   cellular motility, especially the cilium movement of the epithelial
   cells, were also significantly enriched, and most of them ranked in the
   top 20 among all the enriched terms according to their significance
   level (Supplementary Fig. [91]S5-[92]A1, marked with red dots), such as
   the biological process for cilium assembly and movement with GO terms
   of ‘cilium movement’, ‘axoneme assembly’, ‘cilium movement involved in
   cell motility’, ‘outer dynein arm assembly’, ‘epithelial cilium
   movement’, ‘regulation of epithelial cell migration’,
   ‘microtubule-based movement’, the cellular component for the cilium
   with GO terms of ‘cilium’, ‘microtubule’, ‘motile cilium’, ‘axoneme’,
   ‘axonemal central apparatus’, ‘dynein complex’, ‘ciliary basal body’,
   and the molecular function for ‘microtubule motor activity’ and
   ‘microtubule binding’, suggesting exuberant activity of the epithelial
   cells. There was another notable enrichment related to the nervous
   system, including the biological process GO terms ‘ventricular system
   development’, ‘lateral ventricle development’, ‘nervous system
   development’, ‘positive regulation of neuron projection development’,
   ‘brain development’, ‘central nervous system development’,
   ‘sphingolipid metabolic process’ and ‘regulation of neuron apoptotic
   process’ (downregulated), cellular component GO terms ‘synaptic cleft’,
   ‘synapse’ and ‘neuron projection’, and molecular function GO term
   ‘neurotransmitter: sodium symporter activity’, implying the involvement
   of the nervous system in the modulation process after an infestation of
   Polydora. The final overrepresented GO terms were those related to
   energy metabolism and were similar to the other two tissues, including
   glycometabolism and ATP catabolism and anabolism related terms (e.g.,
   ‘nucleoside diphosphate kinase activity’, ‘nucleoside diphosphate
   phosphorylation’, ‘nucleoside triphosphate biosynthetic process’ and
   ‘ATPase activity’). All the relevant GO terms are listed in
   Table [93]1, and those ranked in the top 20 are dotted with different
   colours in Supplementary Fig. [94]5S-[95]A1.

   In the central mantle, there was a similar distribution of the enriched
   GO terms for upregulated genes with terms related to calcium ion
   activity, shell organic matrix, energy metabolism, nervous system and
   cellular motility gathering, as well (Table [96]2; terms ranking in the
   top 20 are dotted with different colours in Supplementary
   Fig. [97]5S-[98]A2). However, differences still existed between the two
   regions. For example, more terms related to the shell organic matrix
   were enriched in the central mantle and some showed upregulation
   (Cellular component: ‘proteinaceous extracellular matrix’,
   ‘extracellular region’ and ‘extracellular exosome’. Biological process:
   ‘protein O-linked glycosylation’ and ‘protein glycosylation’. Molecular
   function: ‘glycoprotein 3-alpha-L-fucosyltransferase activity’,
   ‘alpha-(1- > 3)-fucosyltransferase activity’, ‘fucosyltransferase
   activity’), while some also showed downregulation (Cellular component:
   ‘proteinaceous extracellular matrix’, ‘collagen trimer’ and
   ‘extracellular region’. Biological process: ‘glycosaminoglycan
   metabolic process’). For the terms related to the nervous system, fewer
   showed upregulation (‘neuropeptide signaling pathway’, ‘positive
   regulation of neuron projection development’, ‘neurotransmitter
   transport’ and ‘neuropeptide receptor activity’), but more showed
   downregulation (‘synapse’, ‘glycosphingolipid metabolic process’,
   ‘dendrite morphogenesis’ and ‘neural tube closure’). GO terms related
   to cellular motility were mainly focused on epithelial cell migration
   (‘positive regulation of epithelial cell migration’, ‘colon epithelial
   cell migration’ and ‘regulation of epithelial cell migration’), and
   terms related to cilium movement were not found in the central mantle.
   Different regulation processes probably occurred in the two mantle
   regions suggesting their different tissue functions in shell formation.

Table 2.

   GO terms significantly enriched in the central mantle.
   iD Term p(≤0.05) Type iD Term p(≤0.05) Type
   Upregulated in diseased scallops Related to nervous system
   Related to calcium ion GO:0007218 neuropeptide signaling pathway
   2.18E-04 P
   GO:0051480 cytosolic calcium ion homeostasis 1.28E-03 P GO:0010976
   positive regulation of neuron projection development 6.57E-03 P
   GO:0046686 response to cadmium ion 4.63E-02 P GO:0008188 neuropeptide
   receptor activity 7.47E-03 F
   Related to shell organic matrix GO:0006836 neurotransmitter transport
   8.83E-03 P
   GO:0006493 Related to cell migration
   GO:0018392 glycoprotein 3-alpha-L-fucosyltransferase activity 8.50E-04
   F GO:0010634 positive regulation of epithelial cell migration 2.92E-04
   P
   GO:0046920 alpha-(1->3)-fucosyltransferase activity 1.54E-03 F
   GO:0061580 colon epithelial cell migration 5.24E-04 P
   GO:0008417 fucosyltransferase activity 3.14E-03 F GO:0010632 regulation
   of epithelial cell migration 6.75E-04 P
   GO:0005578 proteinaceous extracellular matrix 1.17E-02 C
   GO:0005576 extracellular region 1.36E-02 C Downregulated in diseased
   scallops
   GO:0006486 protein glycosylation 4.52E-02 P Related to shell organic
   matrix
   GO:0070062 extracellular exosome 4.77E-02 C GO:0005578 proteinaceous
   extracellular matrix 2.66E-03 C
   Related to energy metabolism GO:0030203 glycosaminoglycan metabolic
   process 3.05E-03 P
   GO:0016160 amylase activity 0 F GO:0005581 collagen trimer 2.84E-02 C
   GO:0044245 polysaccharide digestion 0 P GO:0005576 extracellular region
   4.49E-02 C
   GO:0005983 starch catabolic process 1.04E-05 P Related to nervous
   system
   GO:0004396 hexokinase activity 5.07E-05 F GO:0006687 glycosphingolipid
   metabolic process 1.95E-03 P
   GO:0019158 mannokinase activity 5.07E-05 F GO:0048813 dendrite
   morphogenesis 4.54E-03 P
   GO:0004340 glucokinase activity 1.39E-04 F GO:0001843 neural tube
   closure 1.23E-02 P
   GO:0006000 fructose metabolic process 1.39E-04 P GO:0045202 synapse
   1.88E-02 C
   GO:0008865 fructokinase activity 1.39E-04 F Related to immume system
   GO:0005536 glucose binding 2.92E-04 F GO:0034340 response to type I
   interferon 5.87E-04 P
   GO:0001678 cellular glucose homeostasis 3.97E-04 P GO:0032608
   interferon-beta production 9.55E-04 P
   GO:0051156 glucose 6-phosphate metabolic process 3.97E-04 P GO:0006954
   inflammatory response 1.96E-03 P
   GO:0006013 mannose metabolic process 2.15E-03 P GO:0051607 defense
   response to virus 4.21E-03 P
   GO:0006006 glucose metabolic process 1.06E-02 P GO:0005044 scavenger
   receptor activity 5.63E-03 F
   GO:0005975 carbohydrate metabolic process 2.53E-02 P
   [99]Open in a new tab

   Among the downregulated DEGs, many immune-related GO terms were
   detected in all three tissues, especially in hemocytes
   (Tables [100]1–[101]3), and most of them ranked in the top 20 among all
   the enriched terms (Supplementary Fig. [102]5S-B dotted with orange
   colour). For example, four related GO terms (‘inflammatory response’,
   ‘defense response to virus’, ‘innate immune response’ and ‘scavenger
   receptor activity’) were enriched in the edge mantle, and five
   (‘response to type I interferon’, ‘interferon-beta production’,
   ‘inflammatory response’, ‘defense response to virus’ and ‘scavenger
   receptor activity’) in the central mantle, while seven were enriched in
   hemocytes (‘activation of innate immune response’, ‘defense response to
   virus’, ‘positive regulation of defense response to virus by host’,
   ‘innate immune response’, ‘inflammatory response’, ‘immune response’
   and ‘scavenger receptor activity’). The term ‘inflammatory response’
   was also enriched for the three-tissue-shared DEGs (Fig. [103]3A,
   dotted with orange colour). However, no upregulated genes related to
   immunity were found to be significantly enriched, probably indicating
   the declining immune competence in diseased scallops. In addition,
   nervous system-related GO terms were also found enriched in hemocytes,
   with ‘neurotransmitter: sodium symporter activity’ being upregulated
   similarly with the edge and central mantle, while ‘neuropeptide
   signaling pathway’ and ‘neuropeptide receptor activity’ were
   downregulated, indicating the neuromodulation in the immune response of
   M. yessoensis.

Genes related to biomineralization, immune and nervous system

   Many biomineralization-related genes were differentially expressed in
   mantle tissues (Table [104]4), and most of them showed upregulation,
   which suggested a more active shell formation in the infested scallops.
   The functions of these genes were mainly associated with calcium
   binding or transport, such as ‘calbindin’ (CL3204Contig1), ‘calmodulin
   (CL51277Contig1), ‘EF-hand calcium-binding domain-containing protein’
   (CL3228Contig1, CL3161Contig1, CL8277Contig1), ‘sodium/calcium
   exchanger’ (CL36832Contig1) and ‘transient receptor potential cation
   channel’ (CL2149Contig2, CL1979Contig2, CL1618Contig1) and shell matrix
   proteins, such as ‘collagen’ (CL12529Contig1, CL2550Contig1), ‘laminin’
   (CL27Contig2), ‘von Willebrand factor A domain-containing protein’
   (CL41489Contig1, CL721Contig1), ‘sclerostin’ (CL12811Contig1),
   ‘glycoprotein’ (CL55247Contig1,CL4117Contig1), ‘fibronectin
   ‘(CL15874Contig1), ‘glutathione peroxidase’ (CL29217Contig1),
   ‘hephaestin’ (CL12821Contig1), ‘ferritin’ (comp87468_c1_seq2_1),
   ‘tenascin’ (CL27Contig7, CL1236Contig1), ‘chitotriosidase
   ‘(CL50207Contig1), ‘laccase’ (CL6551Contig1), ‘tyrosinase-like
   proteins’ (CL14061Contig1) and ‘carbohydrate sulfotransferase’
   (CL6600Contig1, CL2673Contig1, CL55429Contig1). Most of these genes
   were modulated in only one region of the mantle, implying different
   gene activation in two regions to form the different shell layers.

Table 4.

   Genes related to biomineralization in the mantle tissue.
   Transcript Tissue Fold p(≤0.05) Regulated Description e-value
   Calcium binding or transport
   CL3204Contig1 edge mantle 3.23 1.01E-02 up Calbindin-32 2.00E-83
   CL3204Contig1 central mantle 4.02 1.13E-04 up Calbindin-32 2.00E-83
   CL51277Contig1 edge mantle 5.04 4.41E-02 up Calmodulin 7.00E-34
   CL36832Contig1 edge mantle 4.78 1.55E-02 up Sodium/calcium exchanger 3
   3.00E-168
   CL3228Contig1 edge mantle 5.40 2.47E-02 up EF-hand calcium-binding
   domain-containing protein 1 3.00E-41
   CL3161Contig1 edge mantle 3.13 3.87E-02 up EF-hand calcium-binding
   domain-containing protein 5 2.00E-154
   CL8277Contig1 edge mantle 3.71 2.16E-02 up EF-hand calcium-binding
   domain-containing protein 6 4.00E-175
   CL2149Contig2 edge mantle 3.84 1.56E-02 up Transient receptor potential
   cation channel trpm 3.00E-168
   CL1979Contig2 central mantle 8.45 8.21E-03 up Transient receptor
   potential cation channel subfamily M member 6 8.00E-62
   CL1618Contig1 central mantle 3.95 1.58E-02 up Short transient receptor
   potential channel 3 1.00E-167
   Shell matrix protein
   CL12529Contig1 edge mantle 25.04 6.07E-03 up Collagen alpha-5(VI) chain
   7.00E-60
   CL2550Contig1 central mantle 3.13 1.70E-02 up Collagen alpha-1(XII)
   chain 1.00E-25
   CL27Contig2 edge mantle 4.65 4.53E-03 up Laminin subunit alpha-5
   7.00E-08
   CL721Contig1 edge mantle 3.28 1.49E-02 up Von Willebrand factor A
   domain-containing protein 3A 0
   CL41489Contig1 edge mantle 4.86 1.03E-02 up Von Willebrand factor A
   domain-containing protein 3B 1.00E-18
   CL12821Contig1 central mantle 6.55 1.72E-04 up Hephaestin-like protein
   0
   CL27Contig7 edge mantle Inf 6.25E-03 up Tenascin 8.00E-27
   CL1236Contig1 central mantle 3.26 2.80E-02 up Tenascin-X 1.00E-32
   CL55247Contig1 edge mantle Inf 4.52E-02 up Glycoprotein
   3-alpha-L-fucosyltransferase A 9.00E-45
   CL55247Contig1 central mantle 97.22 3.36E-04 up Glycoprotein
   3-alpha-L-fucosyltransferase A 9.00E-45
   CL4117Contig1 central mantle 9.39 4.80E-03 up
   Glycoprotein-N-acetylgalactosamine 3-beta-galactosyltransferase 1
   7.00E-92
   CL6551Contig1 edge mantle 3.33 3.81E-02 up Laccase-2 8.00E-53
   CL6551Contig1 central mantle 2.61 2.48E-02 up Laccase-2 8.00E-53
   CL14061Contig1 central mantle 2.30 1.79E-02 up Putative tyrosinase-like
   protein tyr-3 3.00E-38
   CL6600Contig1 edge mantle 6.80 2.24E-02 up Carbohydrate
   sulfotransferase 8 5.00E-21
   CL55429Contig1 edge mantle 9.51 2.12E-02 up Carbohydrate
   sulfotransferase 11 2.00E-27
   CL2673Contig1 central mantle −3.53 4.55E-02 down Carbohydrate
   sulfotransferase 9 8.00E-29
   CL12811Contig1 central mantle −3.86 2.44E-03 down Sclerostin
   domain-containing protein 1 1.00E-17
   CL15874Contig1 edge mantle −2.92 2.73E-02 down Fibronectin-like
   1.69E-06
   CL29217Contig1 edge mantle −4.09 4.23E-02 down Glutathione peroxidase
   2.00E-68
   comp87468_c1_seq2_1 central mantle −2.25 2.33E-02 down Soma ferritin
   1.00E-50
   CL50207Contig1 central mantle −13.44 2.09E-02 down Chitotriosidase-1
   3.96E-08
   [105]Open in a new tab

   Many immune-related genes were also modulated in the two mantle regions
   and hemoctyes of the infested scallops (Table [106]5). The functions of
   these genes could be classified into the following categories: I)
   Immune recognition, e.g., ‘lectin’ (CL33366Contig1), ‘toll-like
   receptor (TLR)’ (CL12334Contig1, comp147899_c0_seq5_3), ‘scavenger
   receptor’ (CL63572Contig1), ‘C1q domain containing proteins’
   (CL50918Contig1) and ‘fibrinogen-related proteins’ (CL26062Contig1,
   comp144120_c0_seq2_3). II) Immune effectors, e.g., ‘g-type lysozyme’
   (CL3126Contig2), ‘superoxide dismutase (SOD)’ (comp96557_c0_seq1_1,
   comp119797_c0_seq1_3), ‘catalase’ (CL10534Contig1), ‘big defensin’
   (CL8776Contig1) and ‘heat shock proteins (HSPs)’ (CL1392Contig1,
   CL37312Contig1, CL1984Contig1, CL47124Contig1, CL2340Contig1). III)
   Signal transduction, e.g., ‘tumor necrosis factor (TNF) receptor’
   (CL23037Contig1, comp110989_c0_seq6_2), ‘serine protease’
   (CL61221Contig1) and ‘NF-kappa-B inhibitor’ (comp122698_c0_seq6_2).
   Other genes were included, such as ‘immunoglobulin’ (CL47987Contig1,
   CL19464Contig1), ‘integrin’ (comp127636_c0_seq5_2), ‘interleukin’
   (CL13815Contig1), ‘interleukin receptor’ (CL13604Contig1), ‘stimulator
   of interferon genes protein’ (comp129146_c0_seq2_3) and
   ‘pathogen-related protein’ (comp52531_c0_seq2_1). Unexpectedly, many of
   these genes showed downregulation.

Table 5.

   Genes related to the immune response in the mantle tissue and
   hemocytes.
   Transcript Tissue Fold p(≤0.05) Regulated Description e-value
   Immune recognition
   CL33366Contig1 hemocytes 10.89 7.19E-04 up Lectin 2.00E-50
   CL33366Contig1 edge mantle 20.09 1.63E-03 up Lectin 2.00E-50
   CL33366Contig1 central mantle 8.76 5.93E-09 up Lectin 2.00E-50
   CL26062Contig1 edge mantle 3.58 2.28E-02 up Fibrinogen C
   domain-containing protein 1 1.00E-67
   comp144120_c0_seq2_3 central mantle −5.41 3.75E-02 down Fibrinogen C
   domain-containing protein 1 3.00E-52
   CL12334Contig1 hemocytes −46.78 2.41E-03 down Toll-like receptor 8
   7.00E-19
   CL12334Contig1 edge mantle −36.55 3.47E-03 down Toll-like receptor 8
   7.00E-19
   CL12334Contig1 central mantle −17.80 1.22E-02 down Toll-like receptor 8
   7.00E-19
   comp147899_c0_seq5_3 hemocytes −7.31 2.66E-02 down Toll-like receptor 1
   2.00E-29
   comp147899_c0_seq5_3 edge mantle −8.58 3.53E-05 down Toll-like receptor
   1 2.00E-29
   comp147899_c0_seq5_3 central mantle −5.32 5.59E-05 down Toll-like
   receptor 1 2.00E-29
   CL63572Contig1 hemocytes −3.19 2.13E-02 down Scavenger receptor class F
   member 1 2.00E-13
   CL50918Contig1 hemocytes −10.76 3.59E-02 down Complement C1q tumor
   necrosis factor-related protein 2 4.00E-12
   Immune effectors
   CL3126Contig2 hemocytes 9.86 3.66E-02 up G-type lysozyme 7.44E-138
   comp96557_c0_seq1_1 hemocytes Inf 7.36E-03 up Copper/zinc superoxide
   dismutase 1.33E-07
   comp119797_c0_seq1_3 edge mantle 383.06 3.75E-06 up Copper/zinc
   superoxide dismutase 2.02E-07
   CL10534Contig1 central mantle 2.31 1.68E-02 up Catalase 0
   CL8776Contig1 central mantle 4.04 3.27E-03 up Big defensin 1.00E-35
   CL1392Contig1 hemocytes 9.06 4.25E-02 up Heat shock protein STI1
   1.00E-09
   CL47124Contig1 central mantle 4.07 9.89E-03 up Heat shock 70kDa protein
   12A 4.00E-59
   CL37312Contig1 hemocytes −124.11 1.06E-03 down Heat shock 70kDa protein
   12B 9.00E-92
   CL2340Contig1 central mantle −2.67 1.18E-02 down Heat shock 70kDa
   protein 12B 3.00E-45
   CL1984Contig1 edge mantle −17.45 3.67E-02 down Heat shock protein 90 0
   Signal transduction,
   CL61221Contig1 central mantle 3.04 3.40E-02 up Serine protease 23
   1.00E-09
   comp122698_c0_seq6_2 edge mantle 2.76 3.88E-02 up NF-kappa-B inhibitor
   cactus 1.00E-36
   comp110989_c0_seq6_2 central mantle 10.17 1.38E-02 up Tumor necrosis
   factor receptor superfamily member 11B 8.00E-18
   CL23037Contig1 hemocytes −5.71 1.09E-02 down Tumor necrosis factor
   receptor superfamily member 19L 1.00E-05
   Other related genes
   CL47987Contig1 hemocytes 82.24 2.77E-02 up Immunoglobulin superfamily
   member 10 9.00E-09
   CL47987Contig1 central mantle 65.51 1.08E-03 up Immunoglobulin
   superfamily member 10 9.00E-09
   comp129146_c0_seq2_3 central mantle 5.04 2.10E-02 up Stimulator of
   interferon genes protein 5.00E-15
   comp52531_c0_seq2_1 central mantle 32.68 8.65E-04 up Pathogen-related
   protein 1.00E-33
   CL13604Contig1 edge mantle 3.49 4.51E-02 up Interleukin-6 receptor
   subunit beta 1.00E-06
   CL13815Contig1 edge mantle −7.68 1.13E-02 down Interleukin 17-like
   protein 2.00E-06
   CL19464Contig1 hemocytes Inf 1.84E-02 down Low affinity immunoglobulin
   epsilon Fc receptor 3.00E-21
   comp127636_c0_seq5_2 hemocytes −5.26 4.50E-02 down Integrin beta-6
   3.00E-114
   comp127636_c0_seq5_2 edge mantle −17.01 4.63E-02 down Integrin beta-6
   3.00E-114
   comp127636_c0_seq5_2 central mantle −7.29 2.75E-02 down Integrin beta-6
   3.00E-114
   [107]Open in a new tab

   In addition, some neural genes were differentially expressed in the
   three tissues (Table [108]6). Most of the genes showed upregulation,
   especially for genes for Na^+/Cl^− dependent neurotransmitter
   transporters, such as the ‘GABA transporter’(comp139403_c0_seq2_3),
   ‘glycine transporter’ (CL50297Contig1), ‘aurine
   transporter’(CL31867Contig1), ‘proline transporter’ (CL51047Contig1)
   and ‘creatine transporter’ (CL11023Contig1)in the edge mantle, ‘GABA
   transporter’ (CL6833Contig1) in the central mantle, and ‘glycine
   transporter’ (CL22349Contig1) in hemocytes. However, the gene for
   ‘neurotrypsin’ (CL8026Contig1) was downregulated in all the three
   tissues. Nevertheless, genes for some neurotransmitter receptors had
   different situations. For example, the gene for ‘neuronal acetylcholine
   (ACh) receptor’ (CL31305Contig1) was up modulated in the edge mantle,
   while downregulated in the central mantle (CL42803Contig1) and
   hemocytes (CL35487Contig1). In contrast to the ACh receptor, the gene
   for the ‘5-hydroxytryptamine (5-HT) receptor’ was upregulated in both
   the central mantle and hemocytes but was not detected in the edge
   mantle. Some neuropeptide and neuropeptide receptor genes, such as the
   ‘myomodulin neuropeptides’ (CL32152Contig1) and ‘neuropeptide S
   receptor’ (CL939Contig1), ‘FMRF-amide neuropeptides’ (CL65515Contig1)
   and ‘FMRF-amide receptor’ (comp119939_c1_seq1_2) showed upregulation in
   the central mantle, while ‘FMRF-amide receptor’ (CL24831Contig) was
   downregulated in the hemocytes. There were also other neural-related
   genes were modulated in the three tissues as listed in Table [109]6.
   Genes found in this part (Table [110]6) probably mediated the processes
   of shell repair or immune response through neuromodulation.

Table 6.

   Genes related to the nervous system in the mantle tissue and hemocytes.
   Transcript Tissue Fold p(≤0.05) Regulated Description e-value
   CL11023Contig1 edge mantle 6.62 1.68E-02 up Creatine transporter
   2.00E-117
   CL51047Contig1 edge mantle 3.63 2.72E-02 up Proline-rich transmembrane
   protein 1 3.00E-07
   comp139403_c0_seq2_3 edge mantle 3.08 1.97E-02 up Sodium- and
   chloride-dependent GABA transporter 3 2.00E-75
   CL50297Contig1 edge mantle Inf 2.40E-02 up Sodium- and
   chloride-dependent glycine transporter 1 1.00E-123
   CL31867Contig1 edge mantle 4.80 2.67E-02 up Sodium- and
   chloride-dependent taurine transporter 1.00E-18
   CL6833Contig1 central mantle 2.32 4.72E-02 up Sodium- and
   chloride-dependent GABA transporter 2 2.00E-140
   CL41638Contig1 central mantle 2.67 1.12E-02 up Synaptic vesicular amine
   transporter 2.00E-170
   CL29255Contig1 hemocytes 57.78 1.00E-02 up Sodium-dependent neutral
   amino acid transporter B(0)AT2 4.00E-168
   CL22349Contig1 hemocytes 82.69 2.19E-02 up Sodium- and
   chloride-dependent glycine transporter 2 4.00E-99
   CL32152Contig1 central mantle 4.04 2.65E-02 up Myomodulin neuropeptides
   2 9.00E-09
   CL939Contig1 central mantle 3.83 1.04E-03 up Neuropeptide S receptor
   6.00E-06
   CL65515Contig1 central mantle 2.36 2.48E-02 up FMRF-amide neuropeptides
   2.00E-29
   comp119939_c1_seq1_2 central mantle 3.59 4.65E-04 up FMRFamide receptor
   2.00E-19
   CL24831Contig1 hemocytes −2.67 2.33E-02 down FMRFamide receptor
   6.00E-19
   CL31305Contig1 edge mantle 2.84 3.59E-02 up Neuronal acetylcholine
   receptor subunit alpha-6 4.00E-81
   CL42803Contig1 central mantle −2.40 3.32E-02 down Neuronal
   acetylcholine receptor subunit alpha-6 3.00E-33
   CL35487Contig1 hemocytes −3.64 3.71E-02 down Neuronal acetylcholine
   receptor subunit alpha-10 2.00E-90
   CL32064Contig1 central mantle 2.49 3.80E-02 up 5-hydroxytryptamine
   receptor 1 5.00E-07
   CL1569Contig2 central mantle 2.26 4.96E-02 up 5-hydroxytryptamine
   receptor 2.00E-96
   CL1569Contig2 hemocytes 15.57 3.47E-02 up 5-hydroxytryptamine receptor
   2.00E-96
   CL8026Contig1 edge mantle −11.23 3.77E-02 down Neurotrypsin 2.00E-12
   CL8026Contig1 central mantle −66.03 3.35E-03 down Neurotrypsin 2.00E-12
   CL8026Contig1 hemocytes −2.43 4.80E-02 down Neurotrypsin 2.00E-12
   CL50512Contig1 edge mantle 2.55 3.87E-02 up Advillin 1.00E-82
   CL2411Contig2 edge mantle 38.46 3.02E-02 up Ninjurin-1 2.00E-08
   CL2411Contig1 hemocytes −10.71 4.72E-02 down Ninjurin-1 7.00E-12
   CL3371Contig1 central mantle −2.53 1.54E-02 down Agrin 6.00E-120
   CL6573Contig1 hemocytes 7.51 4.77E-03 up Leukocyte tyrosine kinase
   receptor 2.00E-80
   CL30510Contig1 hemocytes −3.22 3.07E-02 down Orexin receptor type 2
   3.00E-07
   [111]Open in a new tab

Discussion

   The main objective of this study was to identify molecular changes in
   M. yessoensis caused by the infection of Polydora and explore the
   molecular response mechanism of the scallops. The study mainly focused
   on the transcriptional changes in the tissues of edge mantle (producing
   the periostracal and prismatic layers of the shell), central mantle
   (producing the nacreous layer of the shell) and hemocytes (classic
   immune tissue in scallops), which were directly related to the disease.
   The results showed significant alterations associated with the
   processes of biomineralization and immunomodulation in scallops to
   resist infection.

   Changes related to biomineralization mainly occurred in the edge and
   central mantles. Genes related to calcium binding, calcium transporting
   and shell matrix proteins were significantly differentially expressed
   between healthy and diseased scallops (Table [112]4), and the relevant
   gene functions were also significantly enriched in the GO enrichment
   analysis (Tables [113]1 and [114]2). Calcium carbonate (95–99%) and
   shell matrix proteins (1–5%) are the main components of the
   shell^[115]32,[116]33. Genes related to calcium binding and
   transporting regulate calcium concentration in the shell matrix and
   initiate calcium carbonate deposition, while the nucleation, growth and
   spatial orientation of calcium carbonate crystals are regulated by
   shell matrix proteins^[117]34–[118]39. Notably, most of the involved
   genes, such as genes for calbindin, calmodulin, EF-hand calcium-binding
   domain-containing protein, collagen, laminin, tenascin, laccase,
   tyrosinase-like protein, carbohydrate sulfotransferase and others
   showed obvious upregulation in the present study, suggesting increased
   shell secretion activity in infested scallops to repair the damaged
   shells. This hypothesis was also proved by the extremely active
   epithelial cells, which secrete the organic matrix involved in shell
   formation^[119]40. As the results in Tables [120]1, [121]2 and
   Supplementary Fig. [122]S5-[123]A1 show, functions related to the
   activities of mantle epithelial cells were significantly upregulated
   and enriched, especially for those related to cilium motility, probably
   to accelerate the secretion and transport of matrix proteins needed for
   shell repair.

   Differences existed between the two mantle regions with many genes
   regulated in only one region, indicating the different molecular
   functions of the two regions. For instance, genes for the calmodulin,
   sodium/calcium exchanger, EF-hand calcium-binding domain-containing
   protein, laminin, von Willebrand factor A domain-containing protein,
   fibronectin, glutathione peroxidase, etc. were only found to be
   differentially expressed in the edge mantle, while genes for the short
   transient receptor potential channel, hephaestin-like protein,
   tyrosinase-like protein, sclerostin domain-containing protein,
   ferritin, chitotriosidase, etc. were only in the central mantle,
   implying these genes possibly participate in different shell layers
   formation. In addition, some new shell matrix protein genes were first
   found in scallops in this study, i.e., genes for glutathione peroxidise
   and hephaestin, which were once thought to be nonfunctional or lost in
   the mollusk shell and independently arose in Lingula and
   corals^[124]41, but possibly function in scallops shell repair.
   Therefore, this work provides significant candidate genes for the
   mechanism study of shell formation, including genes involved in the
   different shell layers formation, and more research is needed to verify
   how these genes function.

   Although hemocytes have been proved playing an important role in
   biomineralization in oyster due to their involvement in mineral
   transport as well as the production of the extracellular matrix of the
   shells^[125]28–[126]31, the present study didn’t show obvious changes
   related to biomineralization in the hemocytes by gene differential
   expression analysis and GO and KEGG enrichment analysis. However, many
   genes related to immunomodulation were significantly differentially
   expressed in hemocytes and mantle tissue (Table [127]5) and the
   relevant GO terms were also significantly enriched
   (Tables [128]1–[129]3, Supplementary Fig. [130]5SB). The results
   suggest that hemocytes may mainly play an immune function during
   Polydora infection in scallops and the mantle tissue was simultaneously
   involved in this process. As an invertebrate, scallops lack adaptive
   immunity and have evolved a series of sophisticated strategies to
   recognize and eliminate various invaders by employing a set of
   molecules to participate in immune recognition (e.g., lectins^[131]42,
   galectins^[132]43,[133]44 and scavenge receptors^[134]45), signal
   transduction (e.g., serine protease cascades^[135]46–[136]48,
   TLR^[137]49,[138]50 and TNF^[139]51 signalling pathways) and
   incapacitation and elimination of invaders (e.g., lysozymes^[140]52,
   heat shock protein^[141]53,[142]54, and antimicrobial
   peptides^[143]55,[144]56). Immune-related genes detected in the present
   study (Table [145]5) are involved in all the processes mentioned above,
   which reflected an integrated immune response reaction in the scallops
   for the first time. Unexpectedly, a number of these genes showed
   significant downregulation, such as genes for toll-like receptor,
   scavenger receptor, complement C1q tumor necrosis factor-related
   protein, heat shock protein 70 and 90, etc. Furthermore, immune-related
   GO functions were significantly enriched in the downregulated but not
   upregulated genes, such as ‘inflammatory response’, ‘defense response
   to virus’, ‘innate immune response’, ‘scavenger receptor activity’,
   ‘activation of innate immune response’, etc., which indicated declining
   immunity of the infested scallops. Potential trade-offs between
   biomineralization and immunity has been demonstrated in C. gigas and C.
   virginica^[146]29. The expression of biomineralization-related genes
   was higher in C. virginica than in C. gigas, while expression was lower
   for immunity genes. Competition existed between these two processes.
   The contrasting regulation of biomineralization (upregulated on the
   whole) and immunity (downregulated) genes in the present study agreed
   with the above findings.

   In addition, the nervous system likely played an important role in the
   response process. A number of neural genes were significantly
   differentially expressed in the two mantle regions and hemocytes
   (Table [147]6) and relevant GO functions were significantly enriched,
   as well (Tables [148]1–[149]3, Supplementary Fig. [150]5S), suggesting
   the involvement of neuromodulation in the response process. Recently,
   neural-immune regulation in hemocytes was identified in mollusks; this
   regulation was conducted by neurotransmitters via a
   “nervous-hemocyte”-mediated neuroendocrine immunomodulatory axis
   (NIA)-like pathway^[151]57,[152]58. In the present study, GO terms
   associated with neurotransmitter activity were significantly enriched
   in the mantle and hemocytes. Though genes for neurotransmitter were not
   detected among DEGs, genes for several neurotransmitter transporters
   and receptors were significantly differently expressed, such as
   Na^+/Cl^− dependent neurotransmitter transporters, ACh receptor and
   5-HT receptor, which partly implied the existence of neural
   immunomodulation in M. yessoensis. Few studies about neuromodulation in
   shell formation have been reported. The neural-associated DEGs in the
   mantle tissues could not be excluded from the regulation process of
   shell formation in this study, which offers a new perspective on the
   mechanism of shell formation; however, further exploration is needed.

   Finally, scallops likely needed more energy during the response
   process. As we know, glycometabolism is the main pathway to supply
   energy for the organism. The results showed that many
   glycometabolism-related GO terms and pathways, e.g., glucose metabolic,
   fructose metabolic and mannose metabolic, were significantly enriched
   in the upregulated genes among all three tissues (Tables [153]1–[154]3,
   Figs [155]3, [156]4 and Supplementary Fig. [157]5SA), which indicated
   an exuberant energy metabolism process in the diseased scallops. The
   generated energy was probably used for shell repair, immune reactions
   or other activities. For example, the active activity of epithelial
   cells in mantle tissues inevitably consumed a lot of energy.
   Consequently, more energy was used for survival and less energy
   remained for growth.

Conclusions

   To investigate the molecular response mechanisms of M. yessoensis
   infested by Polydora, a comprehensive and niche-targeted cDNA sequences
   database for diseased scallops was constructed and transcriptional
   changes in the edge mantle, central mantle and hemocytes were first
   detected in this study. The gene differential expression and enrichment
   analysis showed that the infestation of Polydora caused expression
   changes for genes or pathways involved in shell formation and immune
   response. Different mantle regions had different molecular functions
   probably responsible for the formation of different shell layers. An
   intact immune response process from immune recognition to pathogen
   elimination was first detected in scallops. Finally, neuromodulation
   and glycometabolism were involved in the infection response. These
   results provide a better understanding of shell formation and innate
   immune response in scallops and provide valuable resources for the
   genetic selection of disease-resistant scallops.

Materials and Methods

Sample collection and RNA extraction

   Two-year-old healthy and infested M. yessoensis (Fig. [158]1) were
   collected from the Dalian Zhangzidao Sea. The two groups of scallops
   were acclimated in the laboratory for one week prior to the
   experiments. Filtered and aerated seawater was maintained at
   approximately 8 °C, which is within the optimum temperature range for
   their growth. Three scallops as biological replicates from each group
   were sampled, and tissues of the edge mantle, central mantle and
   hemocytes that had a direct correlation with the disease were
   immediately placed in liquid nitrogen and stored at −80 °C. Total RNA
   was isolated from each sample using an RNAprep pure tissue kit
   (Tiangent, China). The quantity and quality of total RNA were
   determined using the NanoDrop2000 spectrophotometer (Thermo Scientific,
   Wilmington, DE, USA) and agarose gel electrophoresis.

cDNA library constructionand Illumina sequencing

   mRNA was purified from total RNA by oligo(dT) magnetic beads, and each
   paired-end cDNA library was generated using the TruSeq RNA Sample
   Preparation Kit (Illumina, Inc., USA), following the manufacturer’s
   protocol. The quality of the libraries was assessed using the Agilent
   2100 Bioanalyzer (Agilent Technologies, USA). Finally, a total of 18
   cDNA libraries were sequenced on the Illumina HiSeq2500 genomic
   sequencing platform.

Sequencing data processing, assembly and annotation

   Raw sequencing data was processed to remove reads containing adapters
   or ambiguous “N” nucleotides (length exceeded 35 bp) and low quality
   reads (length of the bases with a quality score less than 20 exceeded
   30%). The high-quality clean reads were combined and used for
   transcriptome assembly by the Trinity software (version:
   trinityrnaseq_r20131110)^[159]59 with default parameters. The generated
   contigs were clustered into unigenes by the TGICL software^[160]60,
   which were used as the reference sequences for subsequent analyses. The
   raw data have been deposited into the NCBI SRA database (accession
   number: SRP150161).

   All the unigene sequences were aligned to the NCBI non-redundant
   protein (Nr) database, Swiss-Prot database, the Eukaryotic Orthologous
   Groups (KOG) protein database, Gene Ontology (GO) database and the
   Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway database by
   BLASTX with an E-value of E ≤ 1e-5 to obtain the best putative
   functional annotations of each unigene.

Gene differential expression analysis

   The expression levels of the unigenes were determined by mapping the
   reads from the samples to reference unigenes using Bowtie 2 software
   ([161]http://bowtie-bio.sourceforge.net)^[162]61. The FPKM (fragments
   per kilobase per million mapped reads) method^[163]62 was used to
   calculate the expression levels of the unigenes for each sample, which
   eliminated the effects of gene length and sequencing depth on the
   calculation of gene expression. Significance tests for expression
   differences in the mantle and hemocytes between healthy and infested
   scallops were conducted with negative binomial distribution
   hypothesis-testing, and the false discovery rate (FDR) method was
   applied to multiple tested hypotheses to correct the significant levels
   (p-values) and eliminate the influence of random fluctuations and
   errors. After calibration, a fold change of 2-fold and p ≤ 0.05 was set
   as the threshold. Enrichment analyses were performed by mapping all
   DEGs (differentially expressed genes) to GO and KEGG databases. GO
   functions or pathways with the statistical significance (p-value) lower
   than 0.05, which were evaluated by the hypergeometric distribution
   test, were considered to be enriched. The DEGs, enriched GO terms and
   KEGG pathways related to shell formation and immunity were highlighted.

Electronic supplementary material

   [164]Supplementary Tables and Figures^ (1MB, pdf)
   [165]Dataset 1^ (1.1MB, xlsx)

Acknowledgements