mRNA Expression and Biomarker Responses in Perch at a Biomonitoring Site in the Baltic Sea – Possible Influence of Natural Brominated Chemicals

Perch ( Perca ﬂuviatilis ) has been used in biological effect monitoring in a program for integrated coastal ﬁsh monitoring at the reference site Kvädöfjärden along the Swedish east coast, which is a site characterized by no or minor local anthropogenic inﬂuences. Using a set of physiological and biochemical endpoints (i.e., biomarkers), clear time trends for “early warning” signs of impaired health were noted in the perch from this site, possibly as a result of increased baseline pollution. The data sets also showed relatively large variations among years. To identify additional temporal variation in biological parameters, global mRNA expression studies using RNA sequencing was performed. Perch collected in 2010 and 2014 were selected, as they showed variations in several biomarkers, such as the activity of the detoxiﬁcation enzyme CYP1A (EROD), the plasma levels of vitellogenin, markers for oxidative stress, white blood cells count and gonad sizes. The RNA sequencing study identiﬁed approximately 4800 genes with a signiﬁcantly difference in mRNA expression levels. A gene ontology enrichment analysis showed that these differentially expressed genes were involved in biological processes such as complement activation, iron ion homeostasis and cholesterol biosynthetic process. In addition, differences in immune system parameters and responses to the exposure of toxic substances have now been veriﬁed in two different biological levels (mRNA and protein) in perch collected in 2010 and 2014. Markedly higher mRNA expression of the membrane transporter (MATE) and the detoxiﬁcation enzyme COMT, together with higher concentrations of bioactive naturally produced brominated compounds, such as brominated indoles and carbazoles, seem to indicate that the perch collected in 2014 had been exposed to macro-and microalga blooming to a higher degree than did perch from 2010. These results and the differential mRNA expression between the 2 years in genes related to immune and oxidative stress parameters suggest that attention must be given to algae blooming when elucidating the well-being of the perch at Kvädöfjärden and other Baltic coastal sites.


INTRODUCTION
The aquatic environment is a final sink for most pollutants.Strategies for environmental monitoring and risk assessment are vital for maintaining the aquatic ecosystem.In particular, biochemical and physiological responses (biomarkers) are valuable tools to provide information and to assess the overall quality of the environment.Methods utilizing biomarkers have been used extensively to investigate exposure, effects and health status in fish (e.g., Larsson et al., 2003;van der Oost et al., 2003, Lehtonen et al., 2014;Asker et al., 2016;Hylland et al., 2017).For this purpose, biomarker measurement of the coastal perch (Perca fluviatilis) have been performed for 30 years in the Swedish National Monitoring Program.These field surveys are run in three reference sites, which started in Kvädöfjärden in 1988, in Holmöarna in 1995, and in Torhamn in 2001.The reason for selecting reference sites for these investigations was to provide information about possible temporal large-scale changes in the coastal ecosystems that may have been caused by changes over time due to different biotic and abiotic factors.In addition, the purpose was to build a database of background biological variables to be used as reference data in studies in polluted coastal areas (Sandström et al., 2005).Therefore, the coastal reference sites are located in areas without any known local or regional point sources and are away from large freshwater inflows.
A continuous decline in the health of fish from those reference areas has been demonstrated by using a combined biomarker and fish health parameter approach (Sandström et al., 2005;Hansson et al., 2006a;Hanson et al., 2009).The decline in fish health in the reference sites is obvious from the successive increase in EROD activity, elevated levels of white blood cell counts, disturbed plasma ion balances, and decreasing condition factors (CFs; a body mass index).Parameters related to fish reproduction also seem to be negatively affected, where the relative gonad size of perch from Kvädöfjärden has been continuously reduced since the start of the investigation in 1988.Based on the evaluation of the first 20 years of biomarker data from perch in Kvädöfjärden, it was suggested that the clear time trends that have been observed for EROD activity and the gonado somatic index (GSI) are related to increasing exposure to environmental contaminants (Hanson et al., 2009).
It is difficult to find a simple explanation for the indicated deterioration of fish health.The current knowledge about the actual chemical pollution at the reference sites and the resulting toxic body burden are only fragmented.At the Kvädöfjärden reference site, which is the focus of the present study, the fish body-burden of most "classic" pollutants, which are measured by the National Monitoring of Contaminants in Biota, such as DDT, HCHs, HCB, PCBs and heavy metals, have shown generally decreasing trends over time (Bignert et al., 2017).The occurrence of naturally occurring brominated compounds, such as brominated diphenylethers and dioxins, show temporal variations but no clear time trends in perch from Kvädöfjärden (Haglund et al., 2010).Other so-called "emerging pollutants, " such as different PFASs, show generally increasing trends (Holmström et al., 2005;Faxneld et al., 2016;Bignert et al., 2017) in biota such as herring and guillemot eggs in the Baltic Sea.
Many of these chemicals can cause adverse effects individually, but it is difficult to link to any individual chemical or part of a complex mixture of chemical compounds to the biomarker time trends.Nevertheless, chemical pollution is usually caused by a complex chemical cocktail that contains dozens of chemicals acting in concert.It has repeatedly been shown that the joint toxicity of chemical mixtures can be substantial, even if all individual compounds are present at only low, individually nontoxic concentrations (Kortenkamp et al., 2009).
In addition to that, the biomarker data from the perch in coastal sites in the Baltic Sea show clear time trends, some of the biomarkers show large variations between years.To investigate and compare variations in the perch biomarkers between two years, 2010 and 2014, which showed large variations, were selected for further studies.For example, the biomarker data from these 2 years showed that the EROD activities, plasma levels of vitellogenin and gonad sizes were higher in 2010 than in 2014 in female perch, whereas catalase and glutathione reductase activities and the plasma content of calcium were markedly lower in 2010 than in 2014.To broaden the biological toolbox for assessing the environmental impact of contaminants and to provide more information about the possible differences between years stored samples of perch liver from the 2010 and 2014 samplings were selected for RNA sequencing.The purpose was to identify differences in mRNA expression patterns that might explain the observed differences between years that might have been caused by different exposure scenarios.In addition, a non-targeted chemical analysis of muscle samples from perch collected in 2010 and 2014 was initiated.The purpose was to provide possible causation for the observed biological differences by analyzing possible differences in the content of brominated dioxins and dioxin-like compounds.Some of these compounds have been found in (and are likely produced by) primary production organisms, i.e., macro-and microalgae in the Baltic Sea (Haglund et al., 2007(Haglund et al., , 2010;;Malmvärn et al., 2008;Unger et al., 2009).

Sampling
Perch (P.fluviatilis) were caught with gill nets (mesh size: 30-33 mm bar length) by a local fishermen at the Kvädöfjärden site on the Swedish east coast (Figure 1).After being caught, fish were carefully released directly from the net and kept for 2 to 3 days to allow stress parameters to go down to base levels before sampling in corves situated at the sampling site (Hansson et al., 2006a and references therein).The sampling took place the last week in September every year between 1988 and 2017.All perch in this study were sexually mature females.At the sampling day, fish were killed by a sharp blow to the head, and blood was collected from the caudal vein with a heparin-prepared syringe.Fresh blood was used for measurement of the haematocrit, hemoglobin content, glucose levels and to produce blood smears for the blood cell count.Thereafter, blood was centrifuged for 90 s at 6,000 g, and the plasma was collected and stored at −80 • C.After measuring the weight and length, the fish was cut open and the bile collected with a syringe.The liver was excised and weighed, and one piece was shock-frozen in liquid nitrogen for measurement of enzyme activities.From year 2010 three extra pieces of liver were shock-frozen in liquid nitrogen in separate cryotubes for analyses of additional biomarkers or other analyses (e.g., RNA sequencing) and then were stored in the large liquid nitrogen containing tanks at the Department of Biological and Environmental Sciences, University of Gothenburg.Fish were weighed after dissection for the somatic weight (carcass weight).The fish carcass was frozen at −20 • C and then sent for storage to Swedish Environmental Specimen Bank at the Natural History Museum in Stockholm.Ethical permission for the samplings was approved by the local animal committee in Gothenburg, Sweden.

Blood Parameters
Blood smears on glass slides were stained using May-Grunwald stain for 5 min, followed by Giemsa stain solution for 18 min.Slides were then rinsed in water and left to dry.Glass slides were analyzed microscopically; approximately 2000 cells were counted per glass slide under magnification (×400).The numbers of immature red blood cells, thrombocytes, lymphocytes, and granulocytes were calculated and presented as a percentage of the total blood cells counted.The total amount of white blood cells (WBC)was calculated as the sum of all thrombocytes, lymphocytes, and granulocytes and was presented as a percentage of the total blood cells.The erythrocyte volume fractions (haematocrit) were estimated using haematocrit capillary tubes followed by centrifugation of the blood using a haematocrit capillary centrifuge for 2 min and a microhaematocrit reader.The hemoglobin and glucose concentrations in blood were measured using a cuvette system from Hemocue, with assayed hemoglobin (HemoTrol; Eurotrol) and glucose (GlucoTrol-AQ; EuroTrol) as quality controls.

Plasma Electrolytes
Levels of the ions Na + , K + , Cl − , and Ca 2+ in the blood plasma were determined with Convergys ISE comfort Electrolyte Analyzer, Cölbe, Germany.

Preparation of Liver Samples
Liver samples were homogenized (glass/Teflon) in four volumes of 0.1 M Na + /K + -phosphate buffer (pH 7.4) containing 0.15 M KCl as previously described (Förlin, 1980).The homogenate was centrifuged at 10,000 g for 20 min, and the supernatant was re-centrifuged at 105,000 g for 60 min.The supernatant (cytosolic fraction) was aliquoted, and the pellet (containing the microsomal enzymes) was re-suspended in one volume of homogenization buffer containing 20% glycerol.All preparation steps were carried out on ice, and the samples were stored at −80 • C until analyzed.

Biochemical Analysis
Ethoxyresorufin O deethylase (EROD) activity was measured in the liver microsomal fraction according to the method described by Förlin et al. (1994) using rhodamine as a standard.The reaction mixture contained a 0.1 M sodium phosphate buffer (pH 8.0), 0.5 mM ethoxyresorufin, and 25 to 50 ml of sample in a final volume of 2 ml.The reaction was started with the addition of 10 ml of 10 mM NADPH.The increase in fluorescence was monitored at 530 nm (excitation) and 585 nm (emission).
The activities of glutathione reductase (GR), glutathione S-transferase (GST) and catalase (Cat) activity was measured in the cytosolic fraction as previously described by Stephensen et al. (2002) and Sturve et al. (2005) and references therein.For GR the reaction mixture contained 0.1 mM DTNB (1-Chloro-2,4dinitrobenzene) and 12 mM EDTA, the reaction was initiated by the addition of 4 mM GSSG (oxidized glutathione) and the absorbance was read at 415 nm.For GST the reaction mixture contained 2 mM CDNB (5,5 -Dithiobis (2-nitrobenzoic acid), 1 mM GSH (glutathione), 0.1 M Na phosphate buffer (pH 7.5).The absorbance was read at 350 nm.For the measurement of the Cat activity the samples were diluted in a 0.08 M KPO 4 -buffer pH 6.5, and 0.08 M H 2 O 2 was added to initiate the reaction.Reactions were measured at 240 nm.

Statistics
Biomarker data from female perch was analyzed with the Mann-Whitney U-test using IBM SPSS statistics 25.Time trends were tested with Spearman's correlation analysis.

RNA Extraction
Liver samples from female perch (5 samples from 2010 and 5 samples from 2014) were homogenized in lysis buffer [RNeasy mini plus kit (Qiagen) using a TissueLyser (Qiagen)] at 25 Hz for 6 min.The total RNA was isolated according to the manufacture's instruction using 50% EtOH.The RNA quality was assessed using TapeStation (Agilent Technologies, United States), and the RIN quality values ranged between 9.4 and 10.0.

Library Preparation and Sequencing
The 10 individual samples were first barcoded, and a single Illumina TruSeq stranded library was then generated and sequenced two lanes of Illumina HiSeq2500 with a 2 × 126 bp setup.This resulted in 12.9-14.5 million reads generated per sample.

Filtering and Trimming of Data
The first ten bases of the reads were removed using fastx_trimmer v0.0.14 from the FASTX-toolkit1 .The adapter and primer sequences were removed using cutadapt v1.3 (Martin, 2011).Bases below a Phred score quality threshold of 20 were then removed with fastq_quality_filter v.0.0.14, which was also from the FASTX-toolkit.An additional four bases in the 5 end were then removed using fastx_trimmer after the result had been analyzed by fastqc v.0.11.4.

Transcript Assembly and Annotation
The reads from the 10 samples were then assembled independently using Trinity v.2.2.0 (Grabherr et al., 2011) with the digital normalization option and taking the strand specificity into account.The ten assemblies resulted in 53770-64140 transcripts, of which 40187-48272 had an assembly score greater than the optimized threshold determined by transrate v1.0.3 (Smith-Unna et al., 2016) (in Supplementary Table S1).To reduce transcript redundancy caused by sequence variability between different individuals, the transcripts with a good assembly score were clustered using CD-HIT (Fu et al., 2012) using a sequence identity cutoff of 97%.In addition, only transcripts represented in at least 2 individuals were kept for further analysis which resulted in 65658 sequences.A BUSCO v3.0.2 (Waterhouse et al., 2017) analysis of this dataset identified 58.8% complete, 11.0% fragmented and 30.2% missing transcripts.Transcripts were annotated using Annocript v.1.1.32and the UniRef90 database (accessed in August 2016) using default settings.

Mapping and Statistical Analysis
The raw reads were mapped to the transcript assembly using bowtie2 v.2.2.2 (Langmead and Salzberg, 2012).Prior to mapping, the bases below a quality threshold of 24 were removed using fastq_quality_filter from the FASTX-toolkit.The number of reads mapped to each transcript was determined using samtools idxstats (Li et al., 2009).
The mRNA expression was analyzed using R version 3.4.23 and the package edgeR v. 3.20.9.Low-abundant transcripts with a total count lower than 4 reads across all samples were filtered out.The mRNA expression was then normalized using TMM normalization.Transcripts that were differentially expressed between the 2014 and 2010 populations were identified by fitting a GLM model implemented in edgeR, and the p-values were corrected by using the false discovery rate (FDR).Functional enrichment using gene ontology terms (GOterms4 ) was carried out among all the significantly differentially expressed genes (FDR < 0.05) between 2010 and 2014, using the web service DAVID5 .The UniProt ID of the best-scoring transcript annotation for each transcript was used as an input to the analysis.

Chemical Analyses of Brominated Dioxin-Like Compounds
Fish muscle (without skin and subcutaneous fat) was used for the chemical analysis.Pooled samples were prepared from the 10 same individuals who are used for RNA sequencing, i.e., 5 samples from 2010 and 5 samples from 2014.In each pool, 5 fish with 20 gram muscle from each fish was used, i.e., 100 gram in total to each pool.The muscle was taken from the middle dorsal muscle layer (TemaNord, 1995).The muscle pieces were sampled from perch retrieved from the Swedish Environmental Specimen Bank, at the Swedish Museum of Natural History.
The fish samples were spiked with a suite of 13Clabeled chlorinated dibenzo-p-dioxin and dibenzofuran internal standards (1 ng each), homogenized with sodium sulfate, extracted with organic solvents, and fractionated according to planarity using activated carbon into three fractions containing, amongst others, (1) the bulk of PCBs, (2) the mono-ortho PCBs, and (3) the non-ortho PCBs and polybrominated and polychlorinated dibenzo-p-dioxins and dibenzofurans and dibenzofurans (Haglund et al., 2007).Fractions 2 and 3 contained the planar dioxin-like compounds.Those fractions were further fractionated according to the polarity on Florisil, and three fractions were collected (Norstrom et al., 1988).The two first contained non-polar contaminants, whilst the two latter contained fat and other semi-polar and polar compounds.
The non-polar fractions from Florisil were screened for brominated compounds using comprehensive two-dimensional gas chromatography (GC × GC; Zoex ZX2, Houston, TX, United States) electron-capture negative ion chemical ionization (ECNI) high-resolution mass spectrometry (HRMS; Agilent, 7250, St. Clara, CA, United States).Extracted ion chromatograms (EICs) of the bromide ions were used to trace brominated compounds.The ECNI and electron ionization (EI) spectra were collected for all brominated compounds.Details on the GC × GC-HRMS analyses are given Table 1.An attempt was made to identify as many brominated compounds as possible using a combination of EI-MS library searches (NIST 17 libary; NIST, Gaithersburg, MA, United States) and manual spectra interpretation.
A semi-quantification was performed using the EI data and the most closely matched internal standard (same degree of halogenation) assuming the same molar response factor, and using the sum of all ions detected of each analyte and internal standard.The result of such a semi-quantification is assumed to be within a factor 2-3 of the correct value.The uncertainty in the ratios of the 2010 and 2014 data should, however, be much less.It should be similar to that of quantitative dioxin determinations, i.e., ca.20% uncertainty (expanded uncertainty, k = 2).

Biomarker Time Trends
The time trend for liver EROD activity in perch at Kvädöfjärden has been reported earlier (Sandström et al., 2005;Hansson et al., 2006a;Hanson et al., 2009).In Figure 2, it can be seen that there is still a significant increasing time trend for the longer period from 1988 to 2017, despite a decreasing trend in EROD activity between the years 2009 and 2014.For GSI, the negative time trend reported before is also still significant for the longer time period, regardless of more stable GSI levels in the last 15 years (Figure 3).Other biomarker time trends in perch include an increase in the liver catalase activity since 2012 (Figure 4).In addition, it has also been indicated that plasma calcium ion content and blood glucose levels show clear increasing time trends in the perch at Kvädöfjärden (Larsson et al., 2016, in Swedish).
TABLE 1 | Settings used for the comprehensive analysis of complex samples using two-dimensional gas chromatography with high-resolution mass spectometry.Comparison Between 2010 and 2014

Morphometric Indices
Fish that were collected in the years 2010 and 2014 were selected, since they displayed large differences in physiological parameters and biomarker levels.In total, 52 fish were included; 27 were sampled in 2010 and 25 were sampled in 2014.The sampled fish from the two different years were of the same age (3.6 and 3.5 years, respectively), but the fish sampled in 2014 were slightly longer and had significantly higher weight resulting in a significantly higher CF compared to the fish sampled in 2010.These fish also had significantly larger livers and liver somatic indices.However, fish sampled in 2014 had lower gonad weights, which led to significantly lower gonadal somatic indices.Morphometric indices are displayed in Table 2.

Blood Parameters
The amount of red blood cells, the haematocrit, was higher in fish sampled in 2014 compared to 2010, but this did not reflect in the hemoglobin content, which remained unchanged.Glucose levels were significantly lower in 2014, while blood lactate levels were higher, even though the difference was not statistically significant (Table 3).
The concentrations of all four plasma ions that were analyzed (Na + , K + , Ca 2+ , and Cl − ) were higher in the fish sampled in 2014 compared to those in 2010, even though the difference was only significant for Ca 2+ an Cl − (Table 3).
The blood cell count revealed that the amount of all cells analyzed, except immature red blood cells (iRBC), were lower in fish sampled in 2014 compared to those sampled in 2010.The lower levels in the amount of the total WBC, lymphocytes, granulocytes and thrombocytes were significant for all cell types.However, the amount of iRBC was significantly higher in 2014, which corresponds to the increase in haematocrit (Table 3).

Biochemical Markers
When comparing fish sampled in 2014 with fish sampled in 2010, the results show lower phase 1 detoxification capacity and oestrogenic responses and higher oxidative stress-related responses.Both EROD activities and the vitellogenin levels were significantly lower in 2014 compared to those in 2010.However, the activities of the antioxidant enzymes that were analyzed (GR, GST, and CAT) were all higher in 2014 compared to 2010, with GR and CAT significantly higher (Table 4).

Transcriptome Analysis for 2010 and 2014 mRNA Expression Studies
Global mRNA expression studies using RNA sequencing were performed to identify additional temporal variations in biological parameters by using perch collected in 2010 and 2014.RNA sequencing analysis identified 65538 genes, of which 4803 had a significant difference in mRNA expression levels (FDR < 0.05).Of these, 2769 genes had a higher mRNA expression level and 2034 genes had a lower mRNA expression level in 2014 compared to 2010.
The genes coding for proteins with biological functions associated with the innate immune system, such as the   toll-like receptor, complement components, serum amyloid A and lysozyme, had a higher mRNA expression level in fish collected in 2014 (Supplementary Table S2).Proteins encoded by genes that are involved in egg formation, such as vitellogenin and zona pellucida, had a higher mRNA expression level in 2010.In addition, several lectins had  different mRNA expression levels in 2010 compared to 2014 (Supplementary Table S2).
The data seem to indicate differently expressed genes associated with oxidative stress, such as PPAR-alpha (peroxisome proliferator-activated receptor alpha) and G6PD (glucose-6phosphate dehydrogenase), and possibly indicate the antioxidant protein peroxiredoxin-6 with higher mRNA levels in 2014 than in 2010.The data also show markedly higher expression of mRNA levels for some membrane pumps, especially MATE1 (multidrug and toxin extrusion protein 1), in the 2014 collected fish (Supplementary Table S2).
The data show higher mRNA expression of genes coding for proteins involved in detoxification, especially of CYP 2, such as CYP 2K1, and phase II proteins, such as UPGT (glucuronosyltransferase), in the fish collected in 2010.In addition, the mRNA expression of genes coding for proteins that are involved in iron homeostasis, such as hepcidin, haptoglobin, and ferritin, showed different patterns in 2010 and 2014.For a ranked list of the top genes according to mRNA expression levels, see the Supplementary Data (Supplementary Table S2).
Principal component analysis (PCA) using all of the sequenced transcripts identified differences in the mRNA expression during the separate time points, as the individual perch collected were clearly divided into two groups (Figure 5).

Gene Ontology Enrichment Analysis
Genes with differential mRNA expressions were analyzed for the overrepresentation of GO terms to identify differences in biological pathways and genes with similar biological functions (Table 5).When using the GO enrichment analysis tool (DAVID) against the human database, several groups of biological processes were found (p < 0.001).In total, 124 regulated biological processes were identified, such as the complement activity, the alternative pathway (GO:0006957), iron ion homeostasis (GO:0055072), and response to estradiol (GO:0032355) and cholesterol biosynthetic process (GO:0006695).

Analyses of Dioxin-Like Brominated Compounds
The non-targeted screening for brominated compounds in perch muscle revealed four compound classes, viz.brominated dibenzo-p-dioxins, carbazoles, indoles and methyl-indoles.The total concentrations (pg/g fresh weight) are given in Table 6 and the individual concentrations (pg/g fresh weight) are given in Supplementary Table S3.
Five mono-, di-, and tri-brominated dioxins were identified in the chemical analyses.Three were found at higher concentrations in 2014, while two were higher in 2010.The total concentration of brominated dioxins was 1.8 times higher in 2010 (Table 6).
For the carbazoles, three dibrominated carbazoles were identified.One was only detected in 2010 and two were found at slightly higher concentrations in 2014.The total  a Functional annotation assessed by GO term enrichment analysis using all significantly differentially expressed genes (FDR < 0.05).The results are shown for the GO term biological process from the hierarchic level 4. concentration of the three dibrominated carbazoles was 14.8 times higher in 2010 perch muscle compared to that in 2014 (Table 6).
The chemical analyses also identified nine brominated indoles and four brominated methylindoles in the perch muscle.For the individual concentrations of indoles, all except one were higher in 2014 than in 2010, and the total concentration of brominated indoles was 6 times higher in 2014 than in 2010 (Table 6).For the four brominated methylindole concentrations, two were higher in 2014 and two were slightly higher in 2010.The total concentration of brominated methylindoles was 2.6 times higher in 2014 compared to that in 2010 (Table 6).Generally, the levels of both brominated indoles and methylindoles were higher in the 2014 perch muscle compared to 2010.

Biomarker Time Trends
It was published that studies of the coastal fish monitoring in the reference area of Kvädöfjärden show that the development of fish health status is not satisfactory (Balk et al., 1996;Sandström et al., 2005;Hansson et al., 2006a;Hanson et al., 2009;Larsson et al., 2011).The previously reported time trends for the GSI and EROD activity in the perch from Kvädöfjärden are still significant for the whole period from 1988 to 2017.For the EROD activity, the marked increase seen from year 2002 to 2009 was followed by a marked decrease between the year 2009 and 2014.After that, the EROD activity had varied between a relatively high to low level.For the GSI, the decreasing trend is still significant, mainly because of large relative gonad sizes during the first ten years of the time series, whereas the data for the last 15 years indicates more stable GSI levels.
The symptoms that were observed in the perch in the reference areas were similar to some of the known effects of certain organic pollutants and to the effects in fish that are exposed to complex mixtures of contaminants, e.g., fish living in areas close to industry effluent point sources or close to other large human activities (Andersson et al., 1988;Larsson et al., 2003;Sturve et al., 2005Sturve et al., , 2014;;Hansson et al., 2006bHansson et al., , 2014)).However, the results of the monitoring of environmental pollutants in the reference areas are not in line with the apparent deterioration of health in coastal fish (Bignert et al., 2017).In Kvädöfjärden, the studies show that most of the organic pollutants and metals that were analyzed in perch were reduced or unchanged (Bignert et al., 2017).However, for some pollutants, e.g., mercury and DDE, in the perch's muscle, increases can be seen for a period of approximately 10 years from 2004 and onward.Increased or unchanged levels of a number of classic environmental pollutants indicate that there was a diffuse supply of these in Kvädöfjärden or a release from sediment.The latter could be supported by major changes in the bottom-dwelling community, especially with an increased occurrence of an invasive organism, the Polychaeta worm Marenzelleria (Michalek and Werner, 2012), whose behavior in sediments is suspected to cause the release of old chemicals from the sediment (Hanson et al., 2016, in Swedish;Hanson et al., manuscript).After 2011, the presence of Marenzelleria decreased again.The period of marked increase of the EROD activity seems to coincide with the occurrence of the digging invasive worm Marenzelleria (Hanson et al., 2016, in Swedish).
Studies show that a number of perfluorinated substances, hydroxylated brominated substances and siloxanes increased in Baltic Sea herring (Faxneld et al., 2014(Faxneld et al., , 2016)).Even in eggs from common guillemot from Stora Karlsö, most perfluorinated substances (PFAS) increased (Bignert et al., 2017).The number of chemicals is increasing substantially in our society, and most of them are not monitored today.Exposure to a mixture of different known and unknown chemicals can, therefore, be an important part of the explanation for the observed health effects of the investigated coastal fish.In addition to the chemical pollution, various environmental factors, such as water temperature, food availability, salinity (Förlin et al., 2014, in Swedish;Hanson et al., 2016 in Swedish) and more, for example, algae blooming, can possibly affect fish biochemistry and physiology over time.

Detoxification Enzymes
Yearly monitoring since 1988 reveals, as expected, variations between years.Since it was of interest to compare the results from fish caught when the EROD activity was high with fish from years when the EROD was low, fish from 2010 to 2014 were selected for comparison and for further analyses.
The EROD activity showed significantly different levels in the female perch between the 2 years, with a 2.3 times higher activity in 2010.In our mRNA expression studies, different CYP1 (A, B, and C) genes could be identified, but no significant differential mRNA expression between the two time points could be seen.However, all identified CYP1 genes showed slightly higher mRNA expression levels in 2010 when compared to 2014, which could well be in line with the higher EROD activity level in 2010.
In contrast to the CYP1 genes, a few other phase I and phase II genes showed differential mRNA expressions between 2010 and 2014.The mRNA expression levels for genes coding for CYP2 proteins were 2.7 to 4.4 times higher in perch from the sampling in 2010 than in 2014.CYP2 proteins are important in the detoxification of a variety of endogenous and exogenous compounds.Only two of the expressed CYP 2 genes were identified to any certainty: CYP2K1 and CYP2J2.In fish species, CYP2K1 has been indicated to be involved in lipid metabolism and aflatoxin metabolism, and CYP2J2 is involved in the metabolism of arachidonic acid (Schlenk et al., 2008;Goldstein et al., 2010).For the phase II detoxification hepatic enzyme GST, the mRNA expression analysis showed three different GST genes that were significantly highly expressed in the 2010 fish.GST is a large family of proteins that catalyze the conjugation of reduced glutathione to a wide range of compounds that result in detoxification (Nebert and Vasiliou, 2004;Schlenk et al., 2008).The activity of this phase II enzyme did not show any significant difference between the 2 years.However, our catalytic GST assay is not specific to different GST proteins and may, therefore, not be able to detect the differences that are indicated for the individual GST mRNA levels.The RNA sequencing also identified genes for other detoxification proteins, including significantly higher mRNA expression of UDP glucuronosyltransferase (UDPGT) in 2010.UDPGT enzymes catalyze the conjugation of glucuronic acid to many hydrophobic compounds, thereby facilitating excretion of toxic compounds through bile and urine (Schlenk et al., 2008;Oda et al., 2015).The differences in mRNA expression of the phase I and II genes that are involved in the detoxification metabolism seem to indicate that fish from 2010 express slightly higher relative enzymatic detoxification capacity than fish from 2014.These differences may indicate different scenarios with respect to exposure to potentially toxic chemicals, but the differences are regarded to be relatively small.

Brominated Dioxin-Like Compounds
In the non-targeted chemical analysis of muscle samples from perch from 2010 and perch from 2014, the focus was to identify brominated dioxin-like compounds.No major differences in the levels of brominated dioxins were found between 2010 and 2014.In addition to brominated dioxins, brominated indoles, brominated methylindoles and brominated carbazoles could be found.Both halogenated carbazoles and indoles have been indicated to be aryl hydrocarbon receptor agonists and, thus, to have the potential to upregulate CYP1A and EROD activity (DeGroot et al., 2015;Mumbo et al., 2015).In the present study, the concentrations of the brominated indoles were generally higher in the 2014 collected fish than in the 2010 fish.This difference is not in line with the difference in EROD activity, with higher activities in 2010 than in 2014.This discrepancy may indicate that the brominated indoles are not causing, or are not the major causation for, the difference in EROD activities between 2010 and 2014.
Three different brominated carbazoles were detected.Two of them were slightly higher in 2014 than in 2010.Many halogenated carbazoles are persistent, express dioxin-like properties (Mumbo et al., 2015) and are possibly more potent inducers of the CYP1A/EROD system than are the indoles.The potency of the detected brominated carbazoles to induce the CYP1A/EROD system in the perch from Kvädöfjärden is not known, and it is not known whether the semi-quantified level found in the perch from Kvädöfjärden is high enough to cause effects in the CYP1A/EROD system.It can, therefore, not be ruled out whether brominated carbazoles play any role in causing differences in EROD activities between the 2 years.
Of the brominated compounds, the brominated dioxins, the brominated indole and methylindole groups are probably mainly of natural origin (Gribble, 1999;Haglund et al., 2010).The sources for these compounds can be macro-or micro-algae.Recently, brominated indoles were tentatively identified in red and brown algae from the Baltic Sea (Björklund, 2018).For the brominated carbazoles, the origins are more unclear with the anthropogenic origin, including sources such as the production of dyes, aluminum industry, rubber, coal combustion and more (Guo et al., 2014;Fromme et al., 2018).However, it has been suggested that carbazoles also have natural origins and seem to dominate in sediment from the great lakes (Barbosa et al., 2014;Guo et al., 2017).Since Kvädöfjärden is located away from large anthropogenic activities, it is likely that the major contribution of the carbazoles are natural sources in that area.In addition, halogenated carbazoles were recently tentatively identified in red and brown algae from the Baltic (Björklund, 2018).
It would be of great interest to perform quantitative measurements of the brominated carbazoles and indoles retrospectively in the perch from Kvädöfjärden.This would be a very important task for the future, as it is likely that historic release in the environment of these groups of halogenated compounds are buried in sediments (Guo et al., 2014(Guo et al., , 2017)).In addition, it is likely that there is an on-going release of these compounds into the environment from current anthropogenic and/or natural activities (Guo et al., 2014(Guo et al., , 2017;;Fromme et al., 2018;Wu et al., 2018).

Membrane Pumps and Detoxification
The RNA sequencing also identified that genes coding for proteins in membrane pumps had a markedly higher mRNA expression in 2014 than in 2010.This was especially the case for the membrane pump that denoted multidrug and toxin extrusion protein (MATE), which was expressed more than one hundred times higher in the 2014 perch compared to the 2010 perch.Genes coding for MATE proteins have been characterized in zebrafish (Loncar et al., 2016).Genes for other drug membrane transporters also had a markedly higher mRNA expression in 2014 fish but not to the same extent as for the MATE gene.The membrane pumps/transporters mediate the excretion of drugs and other organic compounds into bile and urine and thus have an important function in the disposition of toxic compounds (Ferreira et al., 2014).The higher mRNA expression of these membrane transporters could indicate that the fish must handle the transportation of chemicals (endogenous or exogenous) out of the cells (and organism) to a larger degree in 2014 than 2010.
In this context, it is very interesting to the present study that the higher content of brominated compounds in 2014 may be reflected in a higher capacity to transport and excrete these brominated compounds in the 2014 fish than in the 2010 fish.We do not know if the indicated membrane pumps can handle the detected brominated compounds.However, it is of particular interest that the MATE transporter is known to be important as the flavonoid transporter (Lee et al., 2014).Plants are rich in flavonoids.This large group of compounds, many of which have probiotic properties, have been studied extensively for their role in food and health (Moon et al., 2006;Viskupciova et al., 2008).Recently, it has been shown that microalga, including cyanobacteria, contain flavonoids (Goiris et al., 2014;Singh et al., 2017).In the Baltic Sea, the blooming of microalga frequently occurs.Although we have not analyzed any polyphenolic flavonoids in our fish, we speculate that the high mRNA expression of the MATE gene may suggest exposure to flavonoids and thus more profound outbreaks of microalga/cynobacterial blooming close to Kvädöfjärden in 2014 than in 2010.
In addition, the mRNA expression level of the gene coding for the enzyme catechol O-methyltransferase domaincontaining protein 1 (COMT-D1) may indicate higher exposure to flavonoids and thus greater algal blooming in 2014 than in 2010.This gene was indicated to have a ten times higher mRNA expression level in the 2014 fish compared to the 2010 fish.Catechol O-methyltransferase (COMT) is a family of enzymes that degrades catecholamines by methylation (such as dopamine and adrenalin) and catechol oestrogens, as well as various drugs and substances having catechol structure.It is known that flavonoids can also be methylated by COMT enzymes (Zhu et al., 1994;Lemanska et al., 2004).This finding strengthens the contention that the 2014 perch may well have been exposed to microalga blooms that produced flavonoids to a higher degree than did the perch that were sampled in 2010.
Algal blooming frequently occurs in the Baltic Sea.HELCOM provides yearly assessments of cyanobacterial blooms in the Baltic Sea in fact sheets (Öberg, 2017).The fact sheets provide an overall assessment of the blooming and a comparison between years in the entire Baltic Sea.They also provide an annual map to indicate the spatial distribution of the blooms.In the fact sheets for 2010, the overall assessment for 2010 blooms "can be considered to be normal" (Hansson and Öberg, 2010) and for 2014, the bloom "can be considered to be slightly above normal" (Öberg, 2014).In addition, the spatial distribution of the blooms show that in 2010, most of the blooms occurred in the southern part of the Baltic Proper, whereas in 2014, the bloom occurred more in the northern part of the Baltic Proper, which is much closer to the coastal reference site Kvädöfjärden.The information in the fact sheets, therefore, support the assertion that fish in 2014 were likely exposed to algal blooming to a higher degree than the fish in 2010 in the reference area of Kvädöfjärden.

Gonad Development and Vitellogenin
GSI was among the parameters that showed differences between the two time points.In the female perch of the age and size that were used in the present study, the gonads were in a state of rapid growth and development at the time when the fish were netted and sampled in the last week of September.In addition, plasma levels of vitellogenin showed differences between the 2 years.Vitellogenin is produced in the liver and is transported to the gonads, where they are incorporated into the eggs as yolk proteins.The GSI was larger and the plasma content of vitellogenin was higher in the 2010 fish.The same pattern with higher mRNA expression in 2010 was also seen for the genes coding for the protein vitellogenin and for the protein zona pellucida, which were also involved in the egg formation.These differences probably reflect that the 2010 fish were in a slightly more advanced stage in the development of eggs in the gonads (i.e., slightly ahead in growth and development of the gonads).
Whether the detected brominated chemicals may play a role in the difference in gonadal development between the 2 years is unknown.Whether these compounds and flavonoids, many of which have profound bioactive properties, affect the well-being of the perch in Kvädöfjärden would be a very important and interesting research topic in the future.Developmental toxicity by flavonoids in zebrafish, including, for example, the inhibition of aromatase and behavioral disturbances, has been indicated in earlier studies (Brugel et al., 2016).Aromatase converts androgens to oestrogens.Since oestrogens regulate vitellogenin synthesis, the lower plasma vitellogenin levels and lower mRNA expression of vitellogenin in 2014 fish compared to 2010 fish may be the result of aromatase inhibition by flavonoids.In the current state, this is only speculation but is an important hypothesis that is well worth pursuing to discover how or if naturally occurring algae blooms that produce more or less bioactive and toxic compounds may affect reproduction physiology in feral fish.

Oxidative Stress
Other differences between the 2 years seem to indicate that there is oxidative stress in fish in 2014 compared to fish that were sampled in 2010.The enzymatic activities of catalase and glutathione reductase were both higher in 2014.Catalase metabolizes hydrogen peroxide to water and oxygen, while glutathione reductase maintains a correct intracellular ratio of reduced and oxidized glutathione.Glutathione is an antioxidant that is oxidized in antioxidant reactions (Halliwell and Gutteridge, 1999).However, these changes were not reflected in changes in the mRNA expression level; on the contrary, the gene coding for catalase in the perch was slightly higher in the 2010 fish.The gene for another important antioxidant enzyme, peroxiredoxin 6, was indicated to have a markedly higher expression level in 2014 than in 2010.Peroxiredoxin enzymes have the same function as catalase, which is to decompose hydrogen peroxide into water and oxygen.Possibly coupled to the indicated higher mRNA expression in of the peroxiredoxin is the higher mRNA expression of thioredoxin-like protein 1, which was also observed in the 2014 samples.This antioxidant can be related to the function of peroxiredoxin by reducing oxidized peroxiredoxin back to its reduced and active form in the reaction with hydrogen peroxide.A previous study by Sheader et al. (2006) show that both perodixiredoxins and thioredoxins are upregulated by cadmium generated oxidative stress.Taken together, these results suggest that the fish in 2014 had higher levels of intracellular hydrogen peroxide, which possibly derived from pro-oxidant-induced ROS generation.Previous studies show that cyanobacterial blooms lead to effects on the antioxidant defense system in fish and leads to oxidation of macromolecules such as lipids (reviewed in Paskerova et al., 2012).
A higher mRNA expression of the gene coding for the enzyme glucose-6-phosphate dehydrogenase (G6PDH) in the pentose phosphate pathway was also seen in the 2014 fish.G6PDH activity has previously been linked to metal exposure and oxidative stress (Grasset et al., 2016).This enzyme is important in the antioxidant system in organisms because the level of NADPH, which is required to maintain the level of intracellular antioxidant glutathione, is required.The higher mRNA expression of G6PDH and the other above-mentioned mRNA expressions indicate that perch from 2014 show a higher exposure of oxidants than those from 2010.The cause for the apparent higher mRNA expression of genes that are involved in oxidative stress is unknown, but the higher mRNA expression may have strengthened the contention of larger exposure to chemicals that originated from algal blooming in 2014 than 2010.
Several antioxidant genes are regulated by the nuclear transcription factor NRF2, which controls over 200 genes that are involved in the antioxidant system.The NRF2 protein is anchored to the cytoskeleton through KEAP1 and is released and translocated into the nucleus during oxidative stress.Our results show no differences in NRF2 mRNA expression between the 2 years, which suggests mild oxidative stress where the translocation of NRF2 present in the cytosol is enough for antioxidant gene control (Kobayashi and Yamamoto, 2005).

Immunology
In the present study, both the innate (non-specific) and the adaptive (specific) immune systems were differentially regulated.Differences in the white blood cell counts seem to indicate a more activated adaptive immune system in 2010.The higher levels of total WBC, lymphocytes, granulocytes and thrombocytes were significant for all cell types.The higher amount of lymphocytes in perch captured in 2010 also coincides with a higher level of mRNA coding for immunoglobulin chains (Ig heavy chain and Ig kappa chain V, V-I, V-III, and V-IV region, Ig mu chain), which possibly indicates a higher amount of antibody-producing B cells within the lymphocytes population.
A differential mRNA expression of genes coding for proteins that are known to be involved in the innate immune system were noted.For example, serum amyloid A, toll-like receptor 5 lysozyme and most of the complement system components (C1q-like protein, C3, C4-A, C5, C7 C8, C9 and factor H-related protein), had higher mRNA expression levels in fish collected in 2014.All of these proteins are involved in bacterial/viral recognition and elimination, thus taking part in the first line of defense against invading pathogens in vertebrates, including fish (Magnadóttir, 2006;Uribe et al., 2011).
Lectins are a group of proteins with the ability to identify and bind to a diverse set of carbohydrate structures and have multiple functions when expressed in animals, plants and microorganisms.In vertebrates and invertebrates, lectins are involved in cell adhesion and carbohydrate recognition in the circulatory system, as well as within cells.Lectins also have the ability to recognize carbohydrate residues on invading microorganisms, which leads to opsonization, agglutination and the activation of the complement system, thus taking part in the immune defense (Vasta et al., 2011).Several of the lectins that are identified in fish are known to be involved in the innate immune system (Russell and Lumsden, 2005;Vasta et al., 2011).Some fish lectins are also known to be induced after stimuli by different infectious agents (Bayne et al., 2001).Lectins such as collectin-12, the C-type lectin, fucolectin-5, galectin, ladderlectin, galactosespecific lectin nattectin and intelectin are found to be differently expressed in perch, some with a higher mRNA expression in 2010 and others with a higher mRNA expression in 2014.
Genes that are involved in iron homeostasis were differentially expressed in the data set, mostly with a higher mRNA expression in perch collected in 2014.Hepcidin is released upon iron overload (Hu et al., 2007), thereby initiating an increase in iron storage, which was expressed more than 500 times higher in 2014 (logFC at 9.1) than 2010.Other genes with a higher mRNA expression in perch collected in 2014 were haptoglobin and haemopexin binding to free hemoglobin and haeme, respectively (Wassell, 2000), in addition to serotransferrin, mitoferrin and the transferrin receptor protein, which were all involved in iron transport and uptake.Earlier studies in fish, such as turbot, channel catfish and sea bass (Dicentrarchus labrax), has shown that genes coding for proteins involved in iron transport and storage are expressed to a higher extent after bacterial infection (Rodrigues et al., 2006;Peatman et al., 2007;Millan et al., 2011).This would then indicate that iron homeostasis could have a suppressive effect on invading bacteria and hamper their growth and survival due to a low amount of available free iron.
The difference in the amount of WBC and the mRNA expression of genes that are involved in the innate immune system between 2010 and 2014 could be due to different types of bacterial or viral infections in the collected perch.However, abiotic factors such as temperature, oxygen, salinity and pH can have stimulating or suppressive effects on the immune system (Bowden, 2008).Seasonal variation in the respiratory burst of phagocytes and the complement system activity has been noted in the common carp (Cyprinus carpio) (Buchtíková et al., 2011), as well as immunosuppressive effects of steroid hormones in male common carp (Watanuki et al., 2002).Feed additives may also modulate the immune system in fish (Akhter et al., 2015).Chemical pollutants that are found in the environment are also known to affect the immune system (Reynaud and Deschaux, 2006).In yellow perch (Perca flavescens) collected in areas that are polluted with heavy metals, PCBs and PAHs, different lysozyme and ceruloplasmin activity were found compared to that of the reference sites (Dautremepuits et al., 2008).In addition, naturally occurring toxins can have an effect on the immune system.Many secondary metabolites, such as the brominated indoles produced by marine organisms, have multiple and interesting bioactive properties, such as inhibiting inflammation, antibacterial and antiviral properties, and antioxidants (Gribble, 2015).In this context, it is interesting that the indole derivative can combat bacterial infection in marine organisms (Yang et al., 2017).Whether the brominated indoles may have this function in the perch from Kvädöfjärden requires further studies to elucidate.However, it is noted that in the fish from 2014, there was a lower number of WBC active when the levels of brominated indoles are higher than those of the fish from 2010, when the content of brominated indoles were lower.

CONCLUSION
In earlier studies, it was suggested that the time trends in certain biomarker responses in perch from the Kvädöfjärden reference site, studied for more than 30 years, can possibly be the result of changes in the baseline pollution.The RNA sequencing study identified approximately 4800 genes that had different mRNA expression levelswhen comparing perch collected from 2 years, 2010 and 2014, selected as they showed relatively large differences in biomarker levels.The results showed that these differentially expressed genes were involved in biological processes such as complement activation, iron ion homeostasis and cholesterol biosynthetic process.The differences in immune system parameters and responses to the exposure of toxic substances could be verified in two different biological levels (mRNA and protein) in perch collected in 2010 and 2014.The perch collected in 2014 showed higher concentrations of bioactive brominated compounds than did perch from 2010.Macro-and microalgae in the Baltic Sea where algal blooming frequently occurs, likely produce these chemicals.The transcriptomic analyses with the indicated higher mRNA expression of genes that are especially related to membrane transporters and methylation of catechols and the differential mRNA expression between the 2 years in genes related to immune and oxidative stress parameters suggest that attention must be given to algae blooming when elucidating the wellbeing of the perch at Kvädöfjärden and other Baltic coastal sites.The results from the present study therefore suggest that algal blooming and naturally produced chemicals may play a role in exposure to contaminants in the perch from this coastal site.

ETHICS STATEMENT
Ethical permission (Dnr 5.8.18-02260/2018; Idnr 001380) for the samplings was approved by the ethical animal committee in Gothenburg, Sweden.

FIGURE 1 |
FIGURE 1 | Map indicating the sampling area, Kvädöfjärden, at the Swedish east coast.

FIGURE 2 |
FIGURE 2 | EROD activity in the liver (pmol/mg protein x min) in female perch collected in the reference site Kvädöfjärden.The points indicate mean values ± standard error of 20-27 fish in each point, and the straight line represents a significant time trend (P < 0.01).

FIGURE 3 |
FIGURE 3 | Gonado somatic index, GSI (%) in female perch collected in the reference site Kvädöfjärden.The points indicate mean values ± standard error of 20-27 fish in each point, and the straight line represents a significant time trend (P < 0.05).

FIGURE 4 |
FIGURE 4 | Catalase activity (µmol/mg prot x min) in female perch collected in the reference site Kvädöfjärden.The points indicate mean values ± standard error of 20-27 fish in each point, and the straight line represents a significant time trend (P < 0.01).

FIGURE 5 |
FIGURE 5 | Principal component analysis (PCA) performed on the 5 individuals from 2010 and the 5 individuals from 2014 using all 65538 transcripts.Individuals were separated between the different years.

TABLE 2 |
Morphometric indices in female perch (Perca fluviatilis) sampled in Kvädöfjärden in 2010 and 2014 a .

TABLE 3 |
Blood parameters in female perch (Perca fluviatilis) sampled in Kvädöfjärden in 2010 and 2014 a .
a Assessment performed on 27 (2010) and 25 (2014) individuals per site and shown as the mean ± standard error.* indicates significant differences between 2010 and 2014, p < 0.05.b Frequency (%) of WBC (white blood cells), lymphocytes, granulocytes, thrombocytes and iRBC (immature red blood cells) of the total number of blood cells counted.

TABLE 4 |
Activities and levels of biochemical markers in female perch (Perca fluviatilis) sampled in Kvädöfjärden in 2010 and 2014 a .Assessment performed on 27 and 25 individuals per site in 2010 and 2014, respectively, except for VTG, where 15 and 12 individuals were used per site in 2010 and 2014, respectively.Data shown as the mean ± standard error, and * indicates significant differences between 2010 and 2014, p < 0.05. a

TABLE 5 |
Biological processes found in gene ontology (GO) terms enrichment analysis a .

TABLE 6 |
Summation of the concentrations (pg/g fresh weight) of different chemical groups in female perch (Perca fluviatilis) sampled in Kvädöfjärden in 2014 and 2010 and the ratio between the two years.