Estimating salinity stress via hsp70 expression in the invasive round goby (Neogobius melanostomus): implications for further range expansion

Species invasions often occur on coasts and estuaries where abiotic conditions vary, e.g. salinity, temperature, runoff etc. Successful establishment and dispersal of non-indigenous species in many such systems are poorly understood, partially since the species tend to show genetic and ecological plasticity at population level towards many abiotic conditions, including salinity tolerance. Plasticity may be driven by shifting expression of heat shock proteins such as Hsp70, which is widely recognized as indicator of physical stress. In this study, we developed a qPCR assay for expression of the hsp70 gene in the invasive round goby (Neogobius melanostomus) and tested the expression response of fish collected from a brackish environment in the western Baltic Sea to three different salinities, 0, 10 and 30. hsp70 expression was highest in fresh water, indicating higher stress, and lower at brackish (ambient condition for the sampled population) and oceanic salinities, suggestive of low stress response to salinities above the population’s current distribution. The highest stress in fresh water was surprising since populations in fresh water exist, e.g. large European rivers and Laurentian Great Lakes. The results have implications to predictions for the species’ plasticity potential and possible range expansion of the species into other salinity regimes.


Introduction
Biological invasions increasingly are reported; however, the mechanisms facilitating establishment and spread under often novel biotic and abiotic conditions are often unresolved and understudied (e.g. Jaspers et al., 2018). Invasions frequently occur in estuarine areas (Preisler et al., 2009), where salinity gradients or fluctuations are present. Salinity thus may be a primary abiotic factor that the individuals and populations need to be able to tolerate and acclimate to when invading new areas. One hypothesis is that invasive organisms may be able to tolerate, and therefore, their populations may display plasticity in, otherwise physiologically stressful conditions via altering gene expression profiles (e.g. Wellband & Heath, 2017).
The round goby Neogobius melanostomus (Pallas 1814) originated from the Ponto-Caspian region and has been a very successful invader in diverse salinity habitats, such as the Baltic Sea, the Laurentian Great Lakes, as well as in several European rivers (Kornis et al., 2012). Within the brackish Baltic Sea, the species has been detected in all sub-basins including the northernmost, low salinity Bothnian Bay as well as in the higher salinity areas of the western Baltic Sea, including its transition zone with the oceanic North Sea (Azour et al., 2015;Kotta et al., 2016;Puntila et al., 2018). The invasive range of N. melanostomus in the Baltic Sea thus covers highly variable salinity conditions ranging from fresh to almost fully saline waters (Fig. 1). There have been speculations, however, that further range expansions of N. melanostomus into the fully oceanic North Sea might be limited by salinity (Kornis et al., 2012). This was supported by an experiment showing that metabolic and osmoregulatory performance and survival are reduced under high salinities (30, and to some extent at salinity 25) compared to salinities 0-20 (Behrens et al., 2017). It is possible that although adult N. melanostomus often tolerate a wide range of abiotic conditions, reproduction, survival, and early life stage growth of populations may be restricted under high salinity (see Green et al., 2019), leading to decreased invasion potential further north-westward into the fully oceanic North Sea (Behrens et al., 2017). Potential phenotypic plasticity in stress tolerance, however, might aid its populations in survival and establishment in novel, suboptimal environments (Moran & Alexander, 2014), and drive differentiation of the round goby Fig. 1 Salinity gradients (modelled bottom salinity in the area (HELCOM 2020)) in the Baltic Sea and in the study area (red square, enlargement shown to the right). The red arrow points at the sample collection site. The Neogobius melanostomus range is based on data presented in Puntila et al. 2018 scaled for HELCOM sub-basins with EU Water Framework Directive waterbodies populations observed in the area (Björklund & Almqvist, 2010).
Heat shock protein (HSP) expression is often reported to be an indicator of abiotic stress and useful as a biomarker (Deane & Woo, 2011). Many HSPs are expressed in cells constitutively, but when encountering stress and disruption in normal cellular processes, cells can amplify the production of a group of HSP families and thereby maintain cell homeostasis. HSPs perform varying roles in the cell and are classified and named according to their molecular size. The HSP70 complex is necessary for translocation and protein folding and has been widely studied as a biomarker of stress in fishes (Iwama et al., 1999), including salinity stress (Larsen et al., 2008;Tine et al., 2010). Further, a growing number of genome scan studies concurrently point to HSP, and hereunder HSP70, as important candidate genes for local adaptation to salinity (e.g. Nielsen et al., 2009;Tine et al., 2010;Metzger et al., 2016;Kokou et al., 2019). HSP70 is coded by two different genes: a constitutive type (hsc70) that remains unchanged or slightly upregulated when cells are exposed to stress, and an inducible type, hsp70, which is upregulated upon stress. Both genes encode proteins that play key roles in the cell as molecular chaperones (Deane & Woo, 2005).
The present study aimed to test whether experimental salinity stress induces expression of the stress gene hsp70, and whether hsp70 expression varies directly with salinity treatment. To do this, we 1) developed qPCR assays for hsp70 and for a reference gene (EF1) for N. melanostomus and 2) tested the qPCR assays to estimate hsp70 gene expression in a Baltic Sea population N. melanostomus acclimated to three salinity levels, fresh water, reference salinity (10) and oceanic salinity (30).

Materials and methods
Ethical permit 2015-15-0201-00546 from the Danish Animal Ethics Committee covered all experiments reported here. All fish were sacrificed at the end of experiment with an overdose of MS-222. No surgical procedures were performed. The experimental conditions did not appear to severely distress any fish involved in the experiments. Human endpoints were used, i.e. fish that showed signs of suffering or distress (cease of feeding, loss of equilibrium, increase of ventilation frequency or lack of normal movements and social interactions) were euthanized by a sharp blow to the head with a priest, after which they were decapitated.
For the development of the species-specific qPCR assay, two male N. melanostomus (11 cm TL) were collected in Helsinki harbour, Finland (N 60.16°, E 24.96°), on December 27th, 2015. The salinity at the surface water at the sampling location was around 5. The range in salinity between 1971 and 2019 has been 1 to 6 (Syke monitoring database). Tissue samples were extracted from brain, muscle and liver tissue and each sample was divided into 2-3 Eppendorf tubes placed at -80°C until mRNA extraction.
For testing the relationships between salinity stress and hsp70 expression, approximately 120 N. melanostomus individuals were caught with fyke nets in the brackish water estuary (salinity 11, annual range 10-14) Guldborgsund (N 54.7°, E 11.85°) in the western Baltic Sea (Fig. 1). Fish were transported live to DTU Aqua, Denmark, in brackish water from the collection site, and held in freshwater (FW; salinity 0) in 700-L holding tanks receiving filtered, recirculated, and well-aerated water [dissolved oxygen (DO): 90-100% air saturation] at 10°C for 4 weeks until experimental acclimation was initiated. The daylight period was held at 12 h. The transition of the fish to freshwater in the lab upon arrival occurred gradually over 48 h. Water quality was measured once a week, ammonium (NH 4 ? ) and nitrite (NO 2 -) levels with a spectrophotometer (7500 Photometer, Palintest instruments Ltd, Gateshead, UK) and pH with a hand-held pH meter (HQ30d Portable pH Meter, Hach, Loveland, US). The inorganic nitrogen never exceeded 0.05 mg l -1 and 0.01 mg l -1 for [NH 4 ? ] and [NO 2 -], respectively, and pH was stable at 7.9. Fish were fed to satiation three times a week with small commercial high-nutrition fish feed pellets on which they all fed after 1-2 weeks of acclimation. A total of 60 adult fish (all males) of comparable size then were allocated into two separate acclimation tanks containing 200L of water with salinity 10 and 30 fish in each tank. Brackish water (salinity 10) was obtained by mixing filtered oceanic sea water (salinity 30) and nonchlorinated tap water (freshwater) and salinity measured with a digital refractometer (Pocket refractometer PAL-06S from Atago). Thirty fish of comparable size were kept in a separate tank with FW. To maintain adequate water quality (i.e. as described above for the freshwater holding system) in the acclimation tanks, the water in each tank was recirculated through 20-L canister filters with a capacity of 3500 L h-1 and 50% exchanged twice a week. After 7d at salinity 10, the salinity in one of the two tanks was increased gradually by 5 units per week to obtain a final treatment salinity of 30 after four weeks. All fish were subsequently kept at their treatment salinity: FW, reference (salinity 10), and oceanic (salinity 30) for an additional four weeks. On day 55, the experiment was terminated; the fish were measured and weighed; and duplicates of gill and liver samples were rapidly dissected from 10 euthanized fish (given an overdose of MS-222) from each salinity treatment. Samples were snap frozen in liquid nitrogen and subsequently stored at -80°C until RNA extraction.
For all fish, mRNA was extracted using a Qiagen RNeasy mini kit following the manufacturer's instructions. Total RNA quality was assessed using a nanodrop spectrophotometry analysis, and RNA was stored in RNase-free water frozen at -20°C. mRNA then was reverse transcribed to cDNA using the Affinity Script qPCR cDNA synthesis kit following the manufacturer's instructions. Samples were diluted tenfold and stored at -20°C until qPCR analysis.
Hsp70 sequence information was not available for N. melanostomus prior to this study, and we therefore developed a qPCR assay using sequence inference from other fishes. hsp70 sequences for the closest related species available (Gobius species, any gobiid, and any Perciformes species) were searched and retrieved from open access databases NCBI (https:// www.ncbi.nlm.nih.gov) and UniProt (https://www. uniprot.org). As result, hsp70 sequences were included for flounder Platichthys flesus (L. 1758)(AF187726.1) and black goby Gobius niger (L. 1758)(AY762086.1). In order to normalize hsp70 expression, primers were developed for the commonly used house-keeping gene EF1 (elongation factor alpha 1). EF1 was chosen since it has been commonly used in similar studies (e.g. Larsen et al. 2008). Sequences for EF1 included flounder Platichthys flesus (EC378684.1), perch Perca fluviatilis L. 1758 (KC513785.1) and a goby (Pomatoschistus sp.)(EF555084.1). BioEdit 7.1.9 Sequence Alignment Editor (Hall, 2013) and ClustalW Multiple Alignment tool were used to align sequences and identify conserved regions. Three primer pairs were designed for both hsp70 and EF1, prioritising gobiid sequences, when available, using NCBI's Primer Blast tool (https://www.ncbi.nlm.nih.gov/tools/primer-blast/index. cgi?LINK_LOC=BlastHome). Primers (Metabion) were tested with single-tube qPCR using mRNA extracted from the two N. melanostomus collected in Helsinki. qPCR reactions included 5 ll of Bio-Rad iQ SYBR Green super mix, primers (1 ll, each), 3 ll MilliQ water, and 1 ll sample. All qPCR reactions were run with the Bio-Rad CFX384 platform with the following cycling parameters: 95°C for 3 min, 95°C for 10 s, 58°C for 20 s and 72°C for 30 s; for 40 cycles. Melting curve analyses were performed for each reaction to confirm specificity.
Based on the qPCR results, single primer set was selected for hsp70 and EF1, respectively (Table 1) and applied to the 60 samples (three treatments; 10 individuals per treatment; two tissues per individual) from the salinity treatment tests, using the same PCR conditions as above. Optimal annealing temperatures differed between the primer sets (59°C for hsp70 and 62°C for EF1) and assays were therefore run separately, using the cycling parameters otherwise as above. To assess assay efficiency, we applied dilution series based on five dilutions of template DNA pooled for five individuals; resulting in relative concentrations of 1, 0.2, 0.04, 0.008, and 0.0016 and estimated their efficiency following the standard approach (Bustin, 2000). Both runs included two non-template controls (Table 2).
Differences in gene expression between salinity treatments were analysed by comparing the relative expression of the target gene in relation to a reference gene (EF1) and the average expression under the reference treatment (salinity 10), while correcting for assay efficiencies, as detailed in Pfaffl et al. (2002). For reference genes, expression levels may differ among environmental conditions, tissues, and developmental stages. However, the selected gene EF1 is expected to be one of the most stable reference genes across fishes (e.g. Mitter et al., 2009), and we therefore chose to analyse hsp70 expression against just this single reference gene. Statistical differences in hsp70 gene expression ratios estimated for individual fish among treatments and tissues were analysed using one-way ANOVA. Differences with P \ 0.05 were considered significant. The differences in variances between groups (salinity treatments and tissue types) were further analysed with an F-test (Snedecor & Cochran, 1989). Due to three pairwise tests in each group, alpha was corrected using the Bonferroni correction (Dunn, 1958) and results were considered significant at P \ 0.017. All statistical analyses were done in MS Excel.
To confirm that qPCR primers targeted candidate genes hsp70 and EF1, target regions were Sanger sequenced (Sanger & Coulson, 1975). Since this approach works best for fragments above 200 bp, but since both qPCR gene products were below (124 bp and 159 bp), we designed a second set of primers to increase amplicon size (Table 1). Thus, prior to Sanger sequencing, the two target regions were amplified using PCR. Each reaction (volume 20 ll) consisted of 10 ll ThermoFisher 2 9 MasterMix Ò , 2 ll of each primer, 4 ll of MilliQ water and 2 ll of sample template. The reactions were run in Bio-Rad S1000 as follows: 98°C for 1 min, 98°C for 10 s, 59 or 61°C, depending on the primer for 15 s, 72°C for 20 s and 72°C for 5 min. The performance of the PCR was evaluated with 2% Agarose Gel runs. The products then were purified using Agencourt Ò AMPure Ò technique, following the manufacturer's instructions. The concentration of DNA in each purified sample was measured using Qubit Fluorometric Quantification (ThermoFisher). The sequencing reaction (total volume 10 ll) consisted of 0.25 ll of ABI BigDye 1.1, 2 ll of ABI 5x buffer, 1.3 ll of primer, and 5.45-5.65 ll of MIlliQ and 0.8-1.0 ll of PCR product, depending on the DNA concentration. The reactions were run on an ABI 3730 DNA analyser for capillary Sanger sequencing with 96°C for 1 min, 96°C for 30 s, 59 or 61°C, depending the primer for 15 s, 60°C for 4 min, for 30 cycles and then at 72°C for 10 min. Sequences were trimmed by eye and NCBI's Primer Blast tool (https://www.ncbi.nlm.nih.gov/tools/ primer-blast/index.cgi?LINK_LOC=BlastHome) was used.

Results
Inspection of the qPCR melting curves, indicated that qPCR only a single product was produced per assay. Reaction efficiency was 0.96 for hsp70 (standard curve R 2 = 0.99), and slightly lower for EF1 (0.89, R 2 = 0.99). Both gill and liver tissues exhibited higher Hsp70 expression when the goby individuals had been exposed to fresh water compared to the reference salinity (10)(gills: F [1,18 ] = 5.289, P = 0.034; liver: F [1,18] = 18.344, P \ 0.001, Fig. 2). Hsp70 expression showed no statistically significant differences between salinities 10 and 30 (gills: F [1,18] = 0.001, There was a tendency for greater variation among individuals that had been in FW and to some degree at salinity 30 compared to 10, but the variation in expression was significantly different only for the gill tissue between FW and the reference salinity (10)(F 9 = 14.492, P = 0.00024). Sequencing the target genes confirmed that qPCR products represented the expected candidate genes. Sequences were deposited to NCBI (Table 1). Alignment of hsp70 sequences for non-target species used for primer design, round goby, and the selected primers is shown in Fig. 3.

Discussion
N. melanostomus of brackish water origin exhibited a similar low level of stress under oceanic conditions (salinity 30) as in the reference brackish water salinity (10), with the latter being close to the natural conditions at the site of collection (Fig. 1) in the western Baltic Sea. In their native range in the Black Sea, Azov Sea, and the Caspian Sea, they likewise exist in a variety of brackish salinities and along the rivers, into near freshwaters (Chotkowski & Marsden, 1999;Reid & Orlova, 2002). Brackish conditions may be the most optimal to at least some populations of N. melanostomus since their growth rate has also been reported higher in brackish salinities than in fresh water, which may explain why N. melanostomus tend to be larger in the Baltic Sea than in the Great Lakes (Karsiotis et al., 2012). However, our results differ from previous studies investigating salinity tolerance of N. melanostomus originating from fresh water populations in the Great Lakes, where almost 100% mortality in salinities 20-30 was observed (Ellis & Macisaac, 2009;Karsiotis et al., 2012;Behrens et al., 2017). That may be because we tested a population from brackish water, whereas Karsiotis et al. tested a population from freshwater in the Great Lakes (Lake Erie).
Considering that N. melanostomus have established successful populations in low salinity parts of the Baltic Sea (salinity 2-3), several European rivers and other freshwater bodies such as the Great Lakes (e.g. Charlebois et al., 2001;Borcherding et al., 2011;Demchenko & Tkachenko, 2017) it appeared curious that freshwater treatment of this brackish water population led to higher hsp70 expression, i.e. was more stressful than was high oceanic salinity (30). This may be due to your origin population, as there is considerable population structure in the round goby (see Brown & Stepien, 2009). In a recent study, applying acoustic telemetry to study seasonal Asterisks indicate significant differences Fig. 3 Sequence alignments indicating primer sites for hsp70 and the Sanger sequencing used to ascertain target gene in fish species used for primer design migration patterns in N. melanostomus in a brackish water fjord (Karrebaek fjord) adjacent to where the present N. melanostomus were caught (Guldborgsund), just three of the tagged 50 fish migrated into a freshwater stream (River Susaa) with exit in Karrebaek fjord to overwinter. Some of the remaining tagged fish stayed in the fjord during the cold period of the year, but most left the fjord to overwinter in the deeper off-shore, brackish waters (mean salinity 6) (Christoffersen et al., 2019). Although other factors than salinity-related physiological stress most likely also influence habitat choice (e.g. predation pressure, prey availability, wave exposure and/or temperature), these results suggest that N. melanostomus from this area may not have a high affinity for freshwater.
Introductions to the Baltic region presumably have occurred on several occasions (as they have to the Great Lakes (Brown & Stepien, 2009)). Still, the nonindigenous populations in this area assumingly trace to brackish water populations in the northern Black Sea, whereas a riverine (i.e. freshwater) population has been suggested to be the donor source for the Great Lakes N. melanostomus populations (Brown & Stepien, 2009). Furthermore, the species has shown strong adaptation to local conditions (sensu Björklund & Almqvist, 2010), which in the western Baltic Sea could lead to more tolerance towards salinity and in the Great Lakes to more tolerance towards freshwater. To test this, a comparable acclimation protocol with subsequent investigation of hsp70 expression should be performed on round gobies from a variety of salinities, e.g. Great Lake population and native populations in the Black Sea.
Although additional testing is required to confirm generality, our results show that N. melanostomus obtained from a brackish water population in the Baltic Sea acclimated to moderate to high salinities and showed limited stress response in terms of expression of hsp70, pointing to that they may be capable of adapting to higher salinities from which they were collected. Notably, there was pronounced inter-individual variation in expression of hsp70, especially in gill tissue from fish in FW. Although the present study was not designed to test for consistency in individual responses, this suggests that some individuals may maintain unaltered (as compared to the reference salinity of 10) hsp70 expression when subjected to freshwater stress, which may reflect plasticity of salinity tolerance across individuals. In other words, some individuals may be more resilient at extreme low salinities. A recent study has likewise suggested that between-individual differences in blood plasma osmolality of N. melanostomus acclimated to a range of salinities, with some individuals maintaining unperturbed osmolality levels at 25 and 30 PSU as compared to 10 PSU, may reflect variation in salinity tolerance (Behrens et al., 2017).
Our study reports on the expression of a single candidate gene analysed in adult N. melanostomus subjected to four weeks of salinity treatment. Further spread of this invasive species into the saline North Sea will likely also be regulated by a combination of the physiological tolerance and phenotypic plasticity of the likely most vulnerable life stages of the species, i.e. eggs and larvae (Green et al., 2019), by abiotic and biotic interactions (e.g. competition/predation), as well as the species' potential for expanding their range through genetically adapting to higher salinities. Also, there can be new short or long-distance introductions as the species is known to spread with shipping (Kotta et al. 2016). To add further strength to predictions about range expansion, it would be of interest to examine the ability of N. melanostomus to respond to environmental heterogeneity dynamics in peripheral versus in central parts of its distribution (Bronnenhuber et al., 2011). We here provide a tool that can be applied broadly to study N. melanostomus abiotic stress response, e.g. to study variation in gene expression within and among local contingents from different parts of an invasion front at both sides of the Atlantic.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.
Funding Open access funding provided by Finnish Environment Institute (SYKE).