Anthropogenically induced environmental changes in the northeastern Adriatic Sea in the last 1 500 years ( Panzano Bay , Gulf of Trieste ) 2

Anthropogenically induced environmental changes in the northeastern Adriatic Sea in the last 1 500 years (Panzano Bay, Gulf of Trieste) 2 Jelena Vidović a *, Rafał Nawrot a , Ivo Gallmetzer a , Alexandra Haselmair a , Adam Tomašových b , 3 Michael Stachowitsch c , Vlasta Ćosović d and Martin Zuschin a 4 5 a Department of Palaeontology, University of Vienna, Althanstrasse 14, 1090 Vienna, Austria 6 b Earth Science Institute, Slovak Academy of Sciences, Dúbravská cesta 9, 84005 Bratislava, Slovak 7 Republic 8 c Department of Limnology and Bio-Oceanography, Center of Ecology, University of Vienna, 9 Althanstrasse 14, 1090 Vienna, Austria 10 d Department of Geology, Faculty of Science, University of Zagreb, Horvatovac 102a, 10 000 11 Zagreb, Croatia 12 *Corresponding author: vidovic.jelena@gmail.com 13 14 Abstract 15


Introduction
The northern Adriatic Sea is densely urbanized and polluted (Lotze et al., 2006;Cozzi and Giani, 2011), and the areas around the Po River, the Venice Lagoon and in the Gulf of Trieste bear the highest pressure (Solis-Weiss et al., 2007;Raccanelli et al., 2009).Panzano Bay, located in the northeastern part of the Gulf of Trieste, is a shallow and sheltered embayment prone to the accumulation of pollutants, with recent anthropogenic pressure coming from agricultural, maricultural, mining and industrial activities (Horvat et al., 1999).
The impact here started nearly 500 years ago with the onset of mercury mining in the hinterland of the bay (Singh and Turner 2009;Covelli et al., 2012), enhanced in the 20th century with intensifying agriculture and mariculture (Aleffi et al., 2006;Rampazzo et al., 2013;Finch et al., 2014), and Published by Copernicus Publications on behalf of the European Geosciences Union.
has continued until present times with increasing port and industrial activities (thermoelectric plant) of the city of Monfalcone (Notar et al., 2001;Pozo et al., 2009).
Such intensive anthropogenic pressures have prompted a growing scientific effort to estimate the effects of pollution on ecosystem composition here.Most attempts have addressed modifications of the marine habitats that occurred in the 20th century using only geochemical (Horvat et al., 1999;Faganeli et al., 2003;Acquavita et al., 2012) or biological proxies (Solis-Weiss et al., 2007).There is however, a growing tendency towards integrated assessments of its present state (Cibic et al., 2007;Melis and Covelli, 2013;Franzo et al., 2015), but until today there have been no multidisciplinary studies assessing the long-term history of the environmental changes in the northeastern Adriatic and thus capturing its preindustrial, undisturbed state.
Such a historical record requires an integrated geochemical and paleoecological approach.Benthic foraminifera, among the most abundant microorganisms in shallow and marginal marine environments, are often used in paleoecological studies.This is because they are highly sensitive to short-term environmental changes (Schönfeld, 2012) and have a high preservation potential, thus providing an excellent temporal record of ecosystem states over the past hundreds to thousands of years (Yasuhara et al., 2012).
The present multidisciplinary study is designed to provide a high-resolution historical record of environmental changes in Panzano Bay, to obtain information on the state of the ecosystem prior to the onset of the most intensive impact and to evaluate the effects of anthropogenic activities in the bay.We obtained geochemical data and foraminiferal assemblages from a 1.6 m long sediment core containing a centennial-scale record of environmental and anthropogenic changes.The core covers approximately the last 500 years, as indicated by radiocarbon calibrated AAR dating of the mollusk shells.
Taking into account the history of potential anthropogenic stressors in Panzano Bay, we assess the following hypotheses: (1) agricultural and maricultural activities produce upcore increases in the concentrations of organic matter, nutrients and trace elements, (2) mining activities and thermoelectric plants generate a progressive enrichment of mercury and persistent organic pollutants and (3) increased pollutants alter the taxonomic composition of foraminiferal assemblages and cause a decline of species abundance and diversity.
To test these hypotheses, we evaluate the pollution in the bay using geochemical proxies (major, minor and trace elements, nutrients, persistent organic pollutants) and quantify the composition and diversity of foraminiferal assemblages.Finally, we reconstruct the chronology of environmental changes in Panzano Bay over the last 500 years and underline the applicability of our results to disturbed shallow coastal ecosystems elsewhere.

Study area
The Gulf of Trieste is a shallow marine basin in the northernmost part of the Adriatic Sea, occupying an area of about 500 km 2 , with an average water depth of 17 m and a maximum of about 25 m (Fig. 1).Seasonal variations of water temperature range from 8 to 24 • C at the surface and 8 to 20 • C in the bottom layer.The salinity of the water in the gulf is typically marine, ranging between 33 and 38.5 ‰ (Ogorelec et al., 1991).
The water enters the gulf in the southeast and continues to the northwest, following the general anticlockwise circulation pattern of the Adriatic Sea.However, the water circulation in the gulf is mostly controlled by tides (range ∼ 0.5 m), winds (strong northeastern Bora) and seasonal variations of freshwater inflow.The Isonzo/Soča and Timavo rivers are the most significant sources of freshwater to the Gulf, with average inflows of about 100-130 m 3 s −1 each (Ogorelec et al., 1991).
The Gulf of Trieste generally shows mesotrophic to oligotrophic conditions, with episodic eutrophication events, accompanied by summer thermal stratification of the water column (Ogorelec et al., 1991;Horvat et al., 1999;Turk et al., 2007).
The main sediment supply comes from the Isonzo River in the north and from the weathering of the Paleogene flysch deposits outcropping along the southern coast of the gulf.The sediment accumulation rates are approximately 1 mm yr −1 in the central part of the gulf and increase to about 2.5 mm yr −1 towards the mouth of the Isonzo River, located in Panzano Bay (Ogorelec et al., 1991;our unpublished data).Surface sediments in this area are mostly silty clays and clayey silts (Zuschin and Piller, 1994) occupied by a high biomass epifauna (Zuschin et al., 1999).
The Gulf of Trieste is affected by many sources of organic and inorganic pollutants, coming from agricultural and industrial activities in the hinterland, as well as from tourist and maricultural activities along its coasts (Notar et al., 2001;Covelli et al., 2006).Panzano Bay is one of the highly impacted areas, with organic pollution coming from mussel farms located along the eastern part of the Gulf of Trieste (Melaku Canu and Solidoro, 2014) and industrial and port areas of the city of Monfalcone, including a thermoelectric plant and several coal, petroleum and other cargo handling piers (Fig. 1).The Monfalcone thermoelectric plant consists of four thermoelectric generator sets powered by coal and fuel oil that became operative in 1965 (The Monfalcone Thermoelectric Plant, 2016).Finally, there is substantial Hg pollution originating from the Idrija mercury mine in the hinterlands and delivered to the bay through the Isonzo river flow (Horvat et al., 1999;Notar et al., 2001).Idrija, situated 50 km west of Ljubljana (Slovenia), was the second-largest Hg mine in the world, operating for nearly 500 years until its definite closure in 1995 (Faganeli et al., 2003;Covelli et al., 2012).During this period, over 5 million tons of Hg ore was mined and many of the residues were spread around the town and its vicinity (Miklavčič, 1999).Most of the Isonzo riverine input of Hg is in particulate form (1500 kg yr −1 ), followed by dissolved Hg at 8.6 kg yr −1 (Faganeli et al., 2003).Dissolved mercury is biogeochemically reactive and tends to accumulate in certain seafood coming from mariculture, presenting social and economic problems for the local population (Faganeli et al., 2003).

Sampling
Three sediment cores (two for sedimentological and one for foraminiferal analyses), 1.6 m long with a diameter of 9 cm, were acquired using an UWITEC piston corer with hammer action (Gallmetzer et al., 2016) from a research vessel in summer 2013.The drilling station is located in the central part of Panzano Bay (45 • 44 7.32 N; 13 • 36 1.74 E) at a water depth of 12.5 m.The uppermost 20 cm of each core was sliced into 2 cm thick intervals in order to attain highresolution data.The rest of the core was sliced into 5 cm thick samples.For statistical purposes and in order to improve compatibility with the lower part of the core, the uppermost 2 cm thick samples were merged into 4 cm thick intervals (reducing the number of samples from 36 to 31).Sediment samples were used to determine grain size, the content of major, minor and trace elements, nutrients and persistent organic pollutants.Core chronology is based on molluscan shells dated by 14 C calibrated amino acid racemization.

Sediment parameters
The grain size of 36 samples was analyzed using a Sedi-Graph III 5120 Particle Size Analyzer for the small frac-tions (< 63 µm) and dry sieving for fractions from 63 µm to > 1 mm.The sediments were classified according to the Shepard's classification (1954).
To analyze elemental concentrations, each sediment sample was gently squeezed to break down aggregates and screened through a polyethylene sieve to remove particles bigger than 1 mm.A part of the screened sediment was dried in an oven at 105 • C until it reached a constant weight (to measure water content).The dried sediment was ground to powder using an agate mortar and pestle before further analyzing the contents of heavy metals and As.The sample (about 0.4 g d.w.) was digested with 8 mL HNO 3 in a microwave oven (Multiwave 3000, Anton Paar, Austria).The digested material was left to cool at room temperature and then filtered through a 0.45 µm nitrocellulose membrane filter.The filtered digestates were diluted with distilled deionized water to 40 mL in a volumetric flask (USEPA, 1994a).The concentrations of the elements (Al, As, Cd, Cr, Cu, Fe, Mn, Ni, Pb, Zn) were determined by inductively coupled plasma atomic emission spectrometry (ICP-AES) (Optima 2100DV, Perkin Elmer, USA) (USEPA, 1994b).Mercury analyses were carried out using atomic absorption spectrophotometry with cold vapor (Analyst 100, Perkin Elmer, USA) (USEPA, 1976).
The quality acceptance protocols required that one blank sample or one certified reference material (BCR-277r estu-arine sediment, Community Bureau of Reference) were digested and analyzed with each batch of fifteen samples.The blank results indicated that the analytical procedure was free from contamination because the concentrations of all metals were below the respective method detection limits.Mean recovery from the certified material ranged between 84 % (Zn) and 103 % (Hg), except for Al (40 %) because the extraction method was not strong enough to break crystalline aluminosilicates.The analytical precision, determined using five replicates of homogenized samples, was estimated to be better than 10 % for all elements.Calibration for inductively coupled plasma atomic emission spectroscopy (ICP-AES) and atomic absorption spectroscopy (AAS) analyses was achieved with prepared external standards via the standard curve approach.Full calibration was performed after every set of 48 samples.The method detection limit for element analysis was defined as 3 times the standard deviation of 10 blank measurements.
Carbon and nitrogen determination was performed following the method of Hedges and Stern (1984), using an elemental analyzer (CHN 2400, Perkin Elmer, USA).The total concentrations (C and TN) were determined on an aliquot of the sample (about 10 mg of dry sediment); the TOC was determined after treatment of another aliquot of the sample with acid vapors.The inorganic fraction (C) was calculated by the difference between the total and organic concentrations.For the instrument calibration, before each daily series of analyses, three replicates of Acetanilide standard were performed.For the quality acceptance protocols one blank sample in every five samples was analyzed.
To analyze the concentrations of persistent organic pollutants, sediment samples were thoroughly mixed, sieved through a 1 mm mesh to remove any debris, and subsequently air dried in the dark at room temperature for 48 h on hexanerinsed aluminum foil.The dry samples were finely ground in an agate mortar.The extraction was performed using a Microwave Sample Preparation System (Multiwave 3000, Anton Paar Graz, Austria), in accordance with the EPA recommendation (method 3546).Two grams of dried sediments were weighed into lined microwave extraction vessels.Then, a 25 mL 1 : 1 acetone/hexane solvent mixture was added.The vessels were then assembled as instructed by the manufacturer and the extraction was conducted during 15 min at 110 • C and 6-10 bars.At the end of the oven program, vessels were cooled to room temperature and the extracts were filtered and rinsed with the same solvent mixture.
The samples were concentrated in a rotating evaporator (Rotavapor-R Buchi, CH), and the sulphur compounds were removed by soaking the extracts with activated copper powder.Purification and fractionation were performed by eluting extracts through chromatography glass columns packed with Silica gel/Alumina/Florisil (4 + 4 + 1 gr).The first fraction, containing PCBs, was eluted with 25 mL of n hexane, whereas the second fraction, containing the PAHs, was eluted with 30 mL of 8 : 2 n hexane/methylene chloride sol-vent mixture (Fossato et al., 1996(Fossato et al., , 1998)).After concentration with a rotary evaporator, the samples were ready for the instrumental analysis.
The identification of PAHs and PCBs was based on matching retention time, and the quantification was obtained from calibration curves established for each compound by analyzing four external standards.Average determination coefficients R 2 of the calibration curves exceeded 0.99 for both PAH and PCB, and the relative standard deviations of the calibration factors were always less than 20 % (average 10 %).The detection limits were 0.05-0.1 and 0.05 ng g −1 for PAHs and PCBs, respectively.Blanks were run for the entire procedure.Recovery and accuracy were validated with IAEA-417 and IAEA-159 sediment sample certified reference materials.Laboratory methods were also validated by intercalibration activities (IAEA, 2001(IAEA, , 2007(IAEA, , 2012)).
Raw concentrations of Hg, Cr, Pb, As, Cd and PCB were compared to Italian sediment quality guidelines (SQG), following the directive D. L. no.172 of 13 October 2015, whereas PAH and Ni threshold concentrations were taken from directive DM 367/2003.Additionally, raw concentrations were compared to two sediment quality criteria used around the world: effects range low (ERL), representing the threshold level below which effects to benthic organisms rarely occur, and effects range medium (ERM), above which effects are likely to occur (Burton, 2002).Finally, trace elements were normalized to a reference element (Al) in order to compensate for grain size and mineralogical effects on the metal variability in samples (Covelli et al., 2006).

Statistical analyses
Before further statistical treatment, 18 environmental variables (grain size and raw concentrations of nutrients and organic and inorganic pollutants) were checked for normality, log transformed when non-normal distribution was detected, and z-standardized to account for different units and scales.Pearson correlations among environmental variables and principal component analysis (PCA) based on these 18 variables were performed to assess their collinearity and stratigraphic distribution.Only clay content was used in the multivariate analyses because other grain size fractions correlate with the percentage of clay.
The total foraminiferal assemblages were used in all analyses by pooling all mesh size fractions for each sample.Species diversity was measured using species richness, the exponential of Shannon entropy, and Fisher's α.The exponent of the Shannon index (H ) corresponds to the number of equally abundant species that would produce the given value of H (Hill, 1973;Jost, 2006).As all three diversity measures strongly depend on the number of sampled individuals, we rarefied our abundance data down to the size of the smallest sample (240 specimens).This procedure was repeated 1000 times and the mean values of species richness, exp (H ) and α with corresponding 95 % confidence intervals were computed across all iterations.
Species-relative abundance data were square-root transformed before multivariate analyses.Non-metric multidimensional scaling (NMDS) based on Bray-Curtis distances was used to visualize gradients in community composition.Rescaling the NMDS space according to the underlying dissimilarity matrix and rotating it with the principal component analysis (PCA) maximized the compositional variation among samples along the first ordination axis (Oksanen et al., 2015).NMDS axis 1 scores thus correspond to the relative position of samples along the main gradient in species composition.The Pearson correlation was used to measure the association between the environmental variables and NMDS axis 1 scores for the subset of samples with available values of elemental concentrations.
Redundancy analysis (RDA) combined with the forward model selection approach was employed to quantify variation in the multivariate composition of foraminiferal assemblages explained by environmental variables.The effects of environmental variables were first tested in single regressions.Most environmental variables, however, show some degree of collinearity, and the forward model selection approach was thus employed to find a subset of factors that maximizes the explanatory power of the model.At each step of the model building algorithm, an environmental variable with the highest partial R 2 was added while considering the effects of the already selected variables, and the significance of the additional contribution was evaluated through a permutation test (10 000 permutations) (Blanchet et al., 2008).
To identify the timing of the major shifts in assemblage composition, we performed chronological clustering, a type of constrained cluster analysis that takes into account the temporal sequence of samples (Birks, 2012).We used the CONISS algorithm (constrained incremental sum of squares agglomerative clustering) implemented in "chclust" function from the "rioja" package (Juggins, 2015).The number of significantly distinct temporal bins was determined by comparing the amount of variance accounted for by a given number of clusters to a random expectation based on the broken stick model (Bennett, 1996).Clustering was performed on the Bray-Curtis distance matrix based on relative abundance data.All statistical analyses were performed in R 3.2.1 (R Core Team, 2015) using "vegan" (Oksanen et al., 2015) and "rioja" (Juggins, 2015) packages.

Chronological framework
Core chronology is based on the radiocarbon calibrated amino acid racemization dating of the bivalve species Varicorbula gibba.First, 13 shells of V. gibba were selected for 14 C dating and analyzed at the Poznan Radiocarbon Laboratory.Radiocarbon ages were converted to calendar years using Calib7.1 (Stuiver and Reimer, 1993), the Marine13 data (Reimer et al., 2013), and a regional marine reservoir correction ( R) in the northeastern Adriatic equal to −61 years (standard deviation = 50 years) (Siani et al., 2000).The extent of amino acid racemization (AAR) in 329 shells was analyzed at Northern Arizona University using reverse-phase high-pressure liquid chromatography (RP-HPLC) and the procedures of Kaufman and Manley (1998).30 specimens of V. gibba were randomly selected from 11, more or less evenly spaced, 4 or 5 cm thick intervals covering the whole core thickness.The rate of AAR was calibrated based on the 13 shells dated with 14 C and on three live collected individuals with the Bayesian model fitting according to Allen et al. (2013).The time-dependent reaction kinetic model with the initial D / L value estimated from data and lognormal uncertainty showed the best calibration between D / L values of aspartic acid and calendar ages.AAR data in 18 shells did not pass screening criteria, and ages of 311 specimens in total were used for core chronology.Median age of V. gibba in the lowermost (145-150 cm) increment is 1616 AD, but the interquartile range of ages of V. gibba in this increment includes shells that died in the 16th century.Therefore, we used median ages to set the chronology of events, but the core effectively captures the past 500 years of environmental history in Panzano Bay (Fig. 2).Radiocarbon calibrated AAR dating of shells of the bivalve species Varicorbula gibba.The dashed line represent the total range of ages of V. gibba, i.e., bounded by the minimum and maximum age, the solid line represents the inter-quartile range of ages, i.e., it is bounded by the 25th and 75th quantiles, and the grey circle refers to the median age.We note that age range and inter-quartile age range increase downcore and death assemblages in the basal core increments are time-averaged to few centuries, most likely due to bioturbational mixing.Therefore, although median age of V. gibba in the lowermost (145-150 cm) increment is 1616 AD, interquartile range of ages of V. gibba in this increment includes shells that died in the 16th century.

Sediment parameters and geochemistry
The grain size distribution is rather homogeneous throughout the Panzano Bay core, with only a slight increase in the contribution of the > 1 mm fraction in the uppermost part (up to 8.9 %).The sediment in the lowermost part of the core is silty clay (50.4-54.5 % clay).Starting from 135 cm toward the upper section of the core, the amount of clay decreases to 43.5-50 % and the sediment changes into clayey silt (Fig. 3).
Principal component analysis based on raw elemental concentrations illustrates the correlation between elements, with the first two axes explaining 74.8 % of the variance of the data (Fig. 4).This approach distinguishes two major groups of elements with different vertical distribution trends (Fig. 3, Table S2), and three elements (Hg, As, C) that do not fall into this grouping and have distinct position in the ordina-tion space.The first group comprises trace (Cr, Cu, Ni, Cd, Mn) and conservative elements (Fe, Al), characterized by positive mutual correlations (Table S3) and a pronounced decrease in the upper 35 cm.The second group includes organic and inorganic pollutants and nutrients whose raw concentrations are stable (Pb, Zn, PCB) or increase only slightly in the lower part of the core (PAH, TN, TOC, P), but sharply increase in the uppermost 35 cm.Concentrations of PAH increase markedly in the upper 35 cm, although it also shows high values at 75 cm.Normalization to Al reveals two pronounced peaks in the concentrations of the elements from the first group: at 125-130 cm core depth and in the uppermost 10 cm.The latter peak is also visible in normalized Pb and Zn values (Fig. 3).The concentration of Hg sharply increases from 12.74 mg kg −1 at 145-150 cm of the core depth, to 44.7 mg kg −1 at 100-130 cm.The Hg values then decrease upcore to a minimum in the surface sediment (8.22 mg kg −1 ).Concentrations of As vary in the lower core (2.14-9 mg kg −1 ) but gradually decrease to 4.3 mg kg −1 in the surface sediment (0-20 cm).Normalization to Al reveals one concentration peak of As in the upper 10 cm; it coincides with the peak of all other trace elements.Total carbon remains constant throughout the core (8-9.35%), except for the lowermost part (3.7 %).

Trends in foraminiferal assemblages
A total of 69 benthic foraminiferal species were identified in the sediments from Panzano Bay, with raw species richness varying between 29 and 41 species in individual samples (26-36 species after rarefaction to 240 individuals; Table S1 in the Supplement).The highest percentage of individuals belongs to the suborder Rotaliina (63-89 %), followed by Miliolina (8-29 %) and Textulariina (1.5-11 %).Relative abundances of suborders are generally stable throughout the core and vary notably only in the uppermost 20 cm (Fig. 5).Diversity is high throughout the core and increases only in the second half of the 20th century.Values of Fisher α index vary from 7.5 in the lower part to 12 in the uppermost sample; the exponential of Shannon index ranges from 14 to 23 and shows the same vertical trend (Fig. 6a).
Epifaunal/infaunal and infaunal taxa dominate the assemblages, having variable abundances in the lower core (late 17th and 18th centuries) and more stable abundances in its upper part (Fig. 6d).In contrast, the number of infaunal species increases distinctly during the 20th century (Fig. 6e).
Foraminiferal species tolerant of both chemical and organic pollution dominate the assemblages (40-60 %), with maximal abundances in the 18th and 19th centuries and the second half of the 20th century.Species known to tolerate only organic pollution make up 18.5-42 % of the assemblage and have opposite temporal trends than the organic/chemical group, with decreasing trends in the above mentioned time intervals (Fig. 6f, Table S1).NMDS ordination and chronological cluster analysis of the assemblages reveal two main groups of samples, with the major shift in relative species abundances starting around 35 cm (Figs.6b, c and 7).This depth approximately corresponds to the late 19th century, ∼ 1860 AD (Fig. 2).In the NMDS space, samples from the lower and middle part of the core (150-35 cm) are associated with positive scores along the first axis and are tightly grouped, indicating relatively homogeneous faunal composition.These assemblages are characterized by dominance of Valvulineria sp., Nonionella sp., non-keeled elphidiids, Ammonia sp., A. tepida and Haynesina depressula (Fig. 5).The 130-135 cm sample (late 17th century) represents an outlier with unusually low abundance of Ammonia sp. and an increased share of epiphytic species.In contrast, the samples from the upper 35 cm of the core are associated with negative scores along the first axis and are widely distributed in the ordination space.The separation of the uppermost part of the core from the rest of the core in the ordination space suggests a continuous, but relatively strong shift in the assemblage composition at the onset of the 20th century (Fig. 7).Here, the major drop in the abundance of Valvulineria sp. and non-keeled elphidiids is accompanied by a growing share of Miliolinella sp., Triloculina sp., Haynesina sp. and Nonion sp. (Fig. 5).
Within each of the two major groups of samples, further clusters are recognizable, defined by the breaks at 85 and 20 cm (Fig. 6c).The lowermost part of the core (150-85 cm) corresponds to the period from ∼ 1600 to ∼ 1800 AD and has variable foraminiferal distribution trends.The middle part of the core (85-35 cm, ∼ 1800 to ∼ 1860 AD) is characterized by more stable foraminiferal abundances and a pronounced decline of the genus Valvulineria.At 35-20 cm (∼ 1860 to 1950 AD) the diversity of the foraminiferal assemblages starts to increase, as do the abundances of epiphytic species.The uppermost sediment (20-0 cm, 1950 until today) is characterized by a further increase in biodiversity and in the abundance of textulariids (Figs. 5 and 6a).

Relationship between foraminiferal assemblages and geochemical proxies
NMDS axis 1 scores are correlated positively with concentrations of Cu, Ni, Cd, Mn, Fe and Al and negatively with TN and PCB (Table S3).The amount of clay does not correlate with axis 1 scores (Table S3).Total nitrogen content explains the highest proportion of variation in assemblage composition (42.4 %) and is the only explanatory variable included in the RDA analysis following the forward model selection procedure (Table 1 S1. less, other elements that closely (positively or negatively) correlate with total nitrogen content explain a significant amount of variation in single RDA analyses (Table 1), including TOC, organic pollutants (PAH and PCB), and trace elements (Mn, Fe, Ni, Cu, Cd and Zn).The assemblages from the topmost sediment layers (20th century) are clearly separated from the middle core assemblages and from assemblages at the base along RDA axis 1 (Fig. S1).This separation reflects the stratigraphic increase in the content of nitrogen, organic carbon and pollutants and the stratigraphic decline in several trace elements (Fig. 3).

The effects of agricultural and maricultural activities
The agricultural use of pesticides and of organic or inorganic fertilizers releases considerable amounts of pollutants into the environment (Campos, 2003;He et al., 2005;Finch et al., 2014).Pesticides contain pollutant elements such as As, Hg, Cr and Pb (Campos, 2003), while fertilizer contamination includes the discharge of macronutrients (N, P, K) and trace elements, including Co, Cu, Fe, Mn and Zn (Finch et  et al., 2007), leading to elevated concentrations of P, N and TOC in the sediment (Holby and Hall, 1991;Hall et al., 1990Hall et al., , 1992;;Mook et al., 2012).
In the Panzano Bay sediments, trace (Cr, Cu, Ni, Cd, Mn) and conservative elements (Fe, Al) have relatively higher concentrations in the lower part of the core and a pronounced decrease in the upper 35 cm.The only discrepancy of this trend occurred in the late 17th century, when the concentrations of all these elements declined, with simultaneous change in grain size (Fig. 3).Although this elemental fluctuation is in phase with grain size variations, absolute changes in sediment grain size are very minor and it is thus difficult to infer changes in environment based on it.Few elements (Cd, Cr and Pb) sporadically and slightly exceed the limits imposed by the Italian SQG (Fig. 3, Table S2).Only Ni has elevated values throughout the core, even when compared to the standards evaluating the effects of trace elements to benthic organisms (ERL and ERM).Ni, Cd and Cr have a high positive correlation with the major constituents of clay minerals (Al and Fe), the main scavengers of trace elements (Romano et al., 2013).This points to a possible grain size and mineralogical effect on the accumulation of these elements throughout the core because the sediment in Panzano Bay is composed of silt and clay fractions (Fig. 3).
In order to account for such natural processes, to identify background levels and to determine excess trace elements related to anthropogenic contamination, normalized values (trace elements / Al ratios) are usually applied (Covelli et al., 2006;Romano et al., 2013).Normalized concentrations of Cr, Cu, Ni, Cd, Zn and Mn in Panzano Bay are low before the 1950s and, together with As and Pb, increased only in the last 30 years (Fig. 3).Such an increase can reflect the rapid development of technology and the intensification of agricultural activities during the 20th century.
Similar vertical trends have been recorded in the Marano Lagoon, located 20 km west of Panzano Bay (Covelli et al., 2013).The Ni concentrations are almost the same in the two areas, while Pb values are slightly higher in Panzano Bay (starting from 1980 until today).The additional source of Pb here could come from industrial or port activities (see below).
The responses of foraminiferal assemblages to elevated trace element concentrations generally include declining species abundance and diversity as well as altered taxonomic composition because more sensitive species die off and more tolerant taxa prevail (Debenay et al., 2000;Coccioni et al., 2009).Foraminifera can assimilate potentially toxic elements by ingesting contaminated detritus or algae, but also by incorporating these elements during test crystallization, leading to test abnormalities (Le Cadre and Debenay, 2006;Frontalini et al., 2009;Martinez-Colón et al., 2009).In foraminiferal assemblages from Panzano Bay, however, no test abnormalities occurred, indicating that the threshold of elemental concentrations for such an impact was never reached during the last 500 years.

Changes in foraminiferal assemblages before the 20th century
During the period 1600-1900 AD, foraminiferal assemblages in Panzano Bay are characterized by stable diversity indices and a high, variable abundance of stress-tolerant genera and species, including Valvulineria sp., Ammonia sp., A. tepida and non-keeled elphidiids (Figs. 5 and 6f).The genus Ammonia (and especially A. tepida) is usually described as being tolerant of all kinds of stress conditions, including organic and heavy metal pollution (Jorissen, 1988;Coccioni et al., 1997;Armynot du Châtelet et al., 2004;Ferraro et al., 2006;Frontalini and Coccioni, 2008).Non-keeled Elphidium species prefer an infaunal mode of life (Murray, 2006) and can be associated with food enrichment of the sediments (Donnici and Serandrei Barbero, 2002;Vidović et al., 2009Vidović et al., , 2014)).These requirements are similar to the infaunal genus Valvulineria, which is adapted to large seasonal variability of organic matter and periodic hypoxic conditions (Jorissen, 1987; Donnici and Serandrei Barbero, 2002;Piva et al., 2008).Moreover, Valvulineria is considered to be representative of environmental conditions prevailing during the "Little Ice Age" (LIA), which include enhanced rainfall, increased fluvial runoff and increased turbidity (Piva et al., 2008).Interestingly, the distinct peak of Valvulineria in Panzano Bay in the early 19th century (Fig. 5) coincides with the maximal abundances of this genus in sediments from the central and south Adriatic.This peak is attributed to one of the coldest and most humid phases of the LIA, characterized by substantially increased river discharge (Piva et al., 2008).The second peak of Valvulineria in Panzano Bay occurred in the 18th century (Fig. 5), thus showing that humid conditions also prevailed in the bay during this period.To conclude, high abundances of Ammonia, non-keeled Elphidium and Valvulineria during the 17th to 19th century suggest fluctuations of the river runoff and organic matter input in Panzano Bay.
In contrast to fluctuations in abundance of these foraminiferal taxa, vertical changes in the concentration of nutrients in the lower and middle part of the core are mild.However, the spacing of increments that were analyzed for nutrients is larger than dense spacing of increments analyzed for the composition of foraminiferal assemblages.In addition, concentrations of nutrients, grain size distribution and vertical changes in foraminiferal assemblages are likely further affected by vertical homogenization by bioturbation, as evidenced by decadal to centennial time averaging of Varicorbula gibba (Fig. 2), thus making difficult to detect the effects of environmental fluctuations occurring at higher (seasonal or yearly) temporal resolution.
The foraminiferal community from Panzano Bay is also highly correlated with several trace elements (Table S3) that accumulated in fine-grained sediments during this period (as discussed above).Although the assemblages show no effects of elevated trace element concentrations in terms of decline of species abundance or diversity, they remain dominated by taxa tolerant of both chemical and organic pollution (Fig. 6f), as observed also in other foraminiferal assemblage in the northern Adriatic (Jorissen, 1987;Donnici and Serandrei Barbero, 2002).These results imply that the community in Panzano Bay has a long history of adaptation to elevated trace element concentrations.

Changes in foraminiferal assemblages during the 20th century
With the onset of the 20th century, the diversity of foraminiferal assemblages starts to increase (mainly with the increase of infaunal taxa, as reported in Naeher et al., 2012), but this trend becomes pronounced only from 1950 AD onwards (Fig. 6a and e).However, the overall assemblage composition in the 20th century changed markedly, relative to the pre-20th century assemblage composition.The uppermost parts of the core show a very strong and directional change in composition lasting up to the present (Fig. 7), whereas the lower and middle parts of the core were characterized by a relatively constant taxonomic composition, with a much smaller multivariate dispersion in NMDS (Fig. 7).Nutrient concentrations (TN, TOC and P) display the same dynamics of increase during this period (Fig. 3).Although upward increasing concentrations of TOC and TN can be partly related to their recycling dynamic, the corresponding increase in pollutants (PAH and PCB) and other observations of major increase in pollutants and organic enrichments in the Gulf of Trieste (Heath et al., 2006) imply that the nutrient increase also reflects intensifying agricultural and maricultural activities in the Gulf of Trieste during the 20th century.The increase in diversity observed in the uppermost part of the core correlates not only with nutrient enrichment (in accord with observations that early stages of eutrophication can increase species richness, Martinez-Colón et al., 2009) but also with higher concentrations of pollutants, thus rather contrasting with the hypothesis that pollution inevitably decreases species richness.Moreover, TN and P in Panzano Bay sediments are similar to the values measured in sediments beneath adjacent mussel farms (Rampazzo et al., 2013;Franzo et al., 2014).Mussel farming here became an important activity by the middle of the 20th century, reaching peak production in 1990 (Melaku Canu and Solidoro, 2014).Intense mussel biodeposition enriches surface sediments underneath the farms in organic matter, causing anoxic conditions (Rampazzo et al., 2013).Nonetheless, the impact of this farming does not significantly alter the overall coastal marine system (Danovaro et al., 2004;Vidović et al., 2009Vidović et al., , 2014)).Rather, strong winds disperse and resuspend surface organic rich sediments over the broader area of the gulf (Franzo et al., 2014) Besides the increase in diversity, the 20th century is marked by a taxonomic change in foraminiferal assemblages: the abundances of Valvulineria sp., Ammonia sp. and nonkeeled elphidiids decrease, whereas Haynesina sp. and epiphytic genera (Miliolinella sp., Triloculina sp.) become more abundant (Figs. 5 and 6d).Additionally, in the second part of the 20th century herbivorous genus Ammonia also slightly increases its abundances, with only A. tepida not following this trend.Relatively higher abundances of these epiphytic, herbivorous genera during this period suggest the presence of seagrasses or macroalgae meadows near the sampling station (Langer, 1993;Mateu-Vicens et al., 2010).Furthermore, a slight shift in the trophic mode of foraminiferal species in the 20th century (increase of herbivorous taxa) indicates enhanced phytoplankton, probably reflecting higher nutrient levels.The distribution of the genera Miliolinella and Triloculina in the Gulf of Trieste has already been related to their feeding preference for diatoms in addition to organic detritus and bacteria (Hohenegger et al., 1993).
Haynesina, another genus commonly found in the studied sediments, is also herbivorous and known to be tolerant of high concentrations of organic matter (Debenay et al., 2001;Armynot du Châtelet et al., 2004;Murray, 2006;Romano et al., 2008).Higher abundances of Haynesina, together with the increase in overall foraminiferal diversity, may be related to the 20th century nutrient enrichment (Fig. 3) because the representatives of this genus indirectly benefit by feeding on enhanced microalgal biomass (Ward et al., 2003).The fau-nal shift in dominance of Valvulineria to Haynesina, together with higher abundance of epiphytic species, suggests milder seasonal variations of river discharge and enhanced microalgal biomass as a consequence of nutrient enrichment.These conclusions are supported by the RDA analysis, pointing to organic enrichment as a key factor controlling the composition of foraminiferal communities in Panzano Bay (Fig. S1, Table 1).
Finally, the increase in abundance of the suborder Textulariina in this uppermost part of the core may be the result of taphonomic processes: agglutinated taxa are susceptible to postdepositional degradation, and the destruction of their tests explains the downcore reduction of their relative abundances (Diz and Francés, 2009).

Idrija mercury mine
The activity of the Idrija mercury mine is well recorded in Panzano Bay sediments.The Hg concentrations during the last 500 years are high and significantly exceed the limits imposed by the Italian SQG, but also ERL and ERM standards (Fig. 3, Table S2).Interestingly, there are some distinct trends: the concentrations are considerably higher during the 18th century and decrease in the 19th and 20th centuries (Fig. 3), corresponding to the history of the mine: the onset of its significant impact on Panzano Bay occurred in the 18th century, when mining activity sharply increased (Covelli et al., 2012).In the early 19th century, metal recovery from the mine improved, thus releasing less Hg into the river (Covelli et al., 2006).
Foraminiferal assemblages in Panzano Bay remained mostly unaffected by these elevated Hg values throughout the observed period.This implies that speciation of mercury and the bioavailability of its species are more relevant than its total concentration (Martinez-Colón et al., 2009;Acquavita et al., 2012).Most of the Hg enters Panzano Bay in particulate (unreactive) form, with only a small fraction of dissolved Hg (Faganeli et al., 2003).This suggests that the mercury species found here are not accessible to foraminiferal assemblages or, if they are bioavailable, their concentrations do not reach values sufficient to produce toxic effects.

The port of Monfalcone
Panzano Bay is also affected by the industrial and port activities of the city of Monfalcone.Although the first port features were established in the early 19th century, the port as it is known today was designed and built in the 1930s (CPM, 2016).In 1965, a thermoelectric plant powered by coal and fuel oil was opened in the industrial area (The Monfalcone Thermoelectric Plant, 2016).One of the main byproducts of coal and oil combustion are persistent organic pollutants: polycyclic aromatic hydrocarbons and polychlorinated biphenyl, contaminants that potentially form highly carcinogenic and mutagenic derivatives (Notar et al., 2001;Pozo et al., 2009).Moreover, the use of antifouling paints in ports produces trace elements as residues, namely Cu and Zn, but also Cd, Cr, Ni and Pb (Singh and Turner, 2009).
The presence of PAH and PCB in Panzano Bay sediments is probably related to industrial activities in the port.Their concentrations are low throughout the core and start to increase in the 19th century and especially since the middle of the 20th century (Fig. 3), corresponding to the opening of the thermoelectric plant.However not even the highest measured concentrations exceed Italian SQG, or ERL and ERM values.In contrast, the concentration peaks of As, Cr, Cu, Ni, Cd, Zn and Pb in the late 20th century reflect not only agricultural sources (see above) but also intensifying port (antifouling paints) and industrial activities (coal and oil burning).
Certain changes in foraminiferal taxonomic composition correlate with the concentrations of persistent organic pollutants (as detected by RDA in Table 1).These include the decrease of the genus Valvulineria and the increase in the abundance of taxa tolerant of chemical pollution (Fig. 6f), primarily the genus Ammonia (Fig. 5).Nevertheless, as the genus Ammonia is also known to tolerate organic enrichment, a synergistic interaction of both processes (chemical and organic pollution) may have caused such community change.

The chronology of environmental changes in
Panzano Bay over the last 500 years Integration of foraminiferal and geochemical proxies, combined with a robust chronological framework based on extensive radiometric dating of mollusk shells, reveals four major phases in the recent history of Panzano Bay.
During the 17th and 18th centuries, the effects of port activities, as well as of agriculture in the surrounding area, on the composition of foraminiferal assemblages are negligible.In the early 18th century, the release of high amounts of mercury into the environment is related to increasing activity at the Idrija mercury mine (Faganeli et al., 2003).These high inputs, however, did not affect foraminiferal communities because the dominant particulate form of mercury is not bioavailable.Environmental conditions in Panzano Bay during this period were probably characterized by fluctuations in the discharge of the Isonzo River and thus in the amount of organic matter input.The foraminiferal community was therefore composed predominantly of stress-tolerant species adapted to such unstable conditions.
During the 19th century, metal recovery at the Idrija mine improved (Covelli et al., 2006) and less mercury was released into the bay.The onset of maricultural activities here dates back to this period, when bivalve farming was established along the eastern coast of the Gulf of Trieste (Melaku Canu and Solidoro, 2014).This period is marked also by the construction of the port of Monfalcone.The effects of bivalve farming as well as of agricultural and port activities remain negligible during this period.
The first half of the 20th century, however, is marked by rapid technological development.In Panzano Bay, agricultural, maricultural and port activities intensified.The associated slight increase of nutrients caused an increase in foraminiferal diversity and a shift in the trophic mode of the species.
In the second part of the 20th century, the Monfalcone thermoelectric plant, powered by coal and fuel oil, became operative.This slightly increased the concentrations of persistent organic pollutants and caused a minor change in the foraminiferal community composition.The nutrient increase that started in the early 20th century extended to this period.As a consequence, the trend of increasing foraminiferal diversity continues until today.

Conclusions
The chronology of changes in the geochemical composition of sediments and foraminiferal assemblages in shallow and sheltered marine embayments of the northern Adriatic reflects agricultural and industrial development, coastal eutrophication and natural variations.Mercury is a major pollutant in the area, whose concentrations during the last 500 years have significantly exceeded Italian sediment qualwww.biogeosciences.net/13/5965/2016/Biogeosciences, 13, 5965-5981, 2016 ity guidelines, ERL and ERM.Surprisingly, these high concentrations have not affected the ecosystem because the mercury species are not bioavailable to foraminifera.The impact of agricultural, maricultural and industrial activities intensified during the second half of the 20th century and is ongoing.This is reflected in increasing concentrations of trace elements and persistent organic pollutants (PAH, PCB), as well as in progressive nutrient enrichment, as it was presumed within the first two hypotheses.However, mining activity did not produce a progressive enrichment of mercury as anticipated in the second hypothesis, due to the improvement of the methods for the metal recovery.Increased pollutants did not cause a decline of species abundance and diversity as suggested by the third hypothesis, as foraminiferal response to such anthropogenic impacts in Panzano Bay are shaped by their long history of adaptation to elevated trace element concentrations, but also as initial stages of eutrophication can positively affect species richness.Consequently, the shift in community composition during the 20th century reflects a combination of factors, including the recorded increase of pollutants, varying natural conditions, but also a preindustrial predisposition of foraminifera here to tolerate trace elemental pollution.
This combination of factors -and therefore our resultsis clearly applicable to many other shallow coastal areas impacted by human activities, which are largely synchronized on a global scale.Finally, our approach points to the importance of using long-term baseline data for evaluating the environmental and ecological status of present-day marine ecosystems.
The Supplement related to this article is available online at doi:10.5194/bg-13-5965-2016-supplement.

Figure 1 .
Figure 1.Study area and location of sampling site in Panzano Bay.
Figure2.Radiocarbon calibrated AAR dating of shells of the bivalve species Varicorbula gibba.The dashed line represent the total range of ages of V. gibba, i.e., bounded by the minimum and maximum age, the solid line represents the inter-quartile range of ages, i.e., it is bounded by the 25th and 75th quantiles, and the grey circle refers to the median age.We note that age range and inter-quartile age range increase downcore and death assemblages in the basal core increments are time-averaged to few centuries, most likely due to bioturbational mixing.Therefore, although median age of V. gibba in the lowermost (145-150 cm) increment is 1616 AD, interquartile range of ages of V. gibba in this increment includes shells that died in the 16th century.

Figure 3 .
Figure 3. Vertical changes in grain size, major, minor and trace elements, nutrients and persistent organic pollutants.

Figure 5 .
Figure 5. Temporal trends in the relative abundance of foraminiferal suborders (following suprageneric classification of Loeblich and Tappan, 1987) and dominant genera and species (represented by > 2 % of individuals in the pooled data).

Table 1 .
. Results of the forward model selection in the redundancy analysis.Proportion of variance explained in the community data (R 2 ), from permutation tests (10 000 permutations were used) are reported for (a) models with a single explanatory variable and (b) for the effects of a second variable added to the model already including total nitrogen.