Benthic Foraminifera as Pollution Biomarkers: a Morphological Approach

Coastal areas are often intervened by anthropic activities, which increase the contamination of toxic agents such as heavy metals. This causes adverse morphological effects on benthic microorganisms, such as foraminifera. This group is one of the most susceptible to environmental deterioration, so they can be used as pollution biomarkers by identifying shell abnormalities. Therefore, 28 sediment samples from northern Chile were analyzed, calculating the Abnormality Index-FAI and its spatio-temporal distributions in benthic foraminifera, as well as the minimum and maximum abnormality percentages and their relationship with heavy metal concentrations, using a generalized non-linear model and a principal component analysis. The results indicated a proportion of abnormal shells within the ranges described for polluted areas conditions, revealing environmental stress conditions. This reflected a change in the environmental conditions in the most recent sediments of the bay. The highest FAI values were observed to the southwest of the bay, caused by the local current system. The species Bolivina seminuda, Buliminella elegantissima, and Epistominella exigua presented a greater number of deformities, allowing them to be used as contamination biomarkers. A significant correlation was found between Ti, Mn, Ni, Va, and Ba with decreased chamber sizes, wrong coiling, scars, and number of abnormality types. This suggests the effect of the particular geochemical conditions of the area on the heavy metals that cause toxic effects on foraminifera. These analyses are an efficient tool for identifying the effects of environmental stress before they occur in higher organisms, mitigating the environmental impact on marine biodiversity.


Introduction
Coastal areas are subject to constant environmental stress, mainly due to anthropogenic contamination from industrialization, urbanization, tourism, agriculture, and aquaculture activities.This generates large amounts of waste containing heavy metals [1,2] and polycyclic aromatic hydrocarbons (PAHs) [3], which generate lethal and non-lethal effects on marine organisms [4] and are one of the main causes of marine biodiversity loss [5].Heavy metals, specifically, show a certain affinity with unicellular organisms, which absorb these agents through the cell membrane through diffusion processes, which control the rate of uptake and binding of proteins within the cell [6].Additionally, in high concentrations, heavy metals generate enzymatic inhibition, causing mortality in sensitive organisms [7].This is why using bioindicators is of the utmost importance, as they allow evaluating the effects of contamination in a given ecosystem, covering not only the presence or absence of certain groups of individuals but also any biotic response produced by environmental stress, evidenced at the molecular, biochemical, physiological, population, community, or ecosystem levels [8].Thus, these biotic responses allow identifying environmental changes that cannot be directly observed [9], generating a tool for monitoring the health status of coastal ecosystems.This allows integrating environmental stress events reflected by qualitative and quantitative changes in the community structure of a certain group, and changes in functional traits within the ecosystem.Among the benthic organisms used as contamination bioindicators are nematodes, mollusks, ostracods, and foraminifera, which are useful biological models to identify environmental changes produced by heavy metals [10,11].
Foraminifera in particular make up a large part of the marine biota and have a cosmopolitan distribution [12,13].They allow for a more detailed and complete characterization of environmental conditions due to the excellent preservation of their calcareous shells over time [14] and the specificity of their microhabitats.This causes a rapid response to environmental changes, reflected by fluctuations in the diversity and abundance of these microorganisms, as well as morphological effects on their shells due to environmental stress conditions.They are, therefore, considered one of the groups with the highest sensitivity to coastal environment deterioration [15,16].Thus, the rapid response of foraminifera, either at the community or morphological level due to adverse environmental conditions, occurs prior to the evidence of negative effects in more complex organisms (invertebrates and vertebrates), allowing their use as early biomarkers of environmental stress.
In Chile, 60 km from Antofagasta is Mejillones Bay, characterized for its low oxygen concentrations associated with the oxygen minimum zone (OMZ) [17].The circulation of the bay is influenced by a source of upwelling in the southern area outside the bay, which induces the flow to move northward and modulates the filling and emptying of the bay [18].This is one of the most productive upwelling centers in the area, generating a high diversity "hotspot" [19][20][21].
At the beginning of the twentieth century, industrial activities increased in the area, most of them focused on aquaculture, fisheries, urbanization, construction of port terminals for the export of minerals, and thermoelectric plants.Despite being a heavily intervened area, few studies have been carried out on invertebrate organisms as bioindicators of adverse conditions in the bay.Among these, we highlight the work of [22], referring to the effect of dredging on the diversity of macrobenthic communities due to the construction of Angamos Port, and an ex situ research carried out by Guiñez et al. [23], who studied the bioaccumulation of heavy metals in the crustacean Emerita analoga in the intertidal zone of the bay.
Therefore, our work is focused on implementing foraminifera as early biomarkers of heavy metal contamination, based on their morphological responses.We aim to identify which are the heavy metals that produce abnormalities in the shells of foraminifera and what types of anomalies they are related with.In addition, we also expect to know what is the spatial and temporal distribution of abnormalities in the Bay of Mejillones and likewise, to be able to compare them with other polluted coastal areas.This contribution to the knowledge regarding the environmental quality of this minimum oxygen zone throughout time and space evidences the relevance of studying these organisms as a useful tool for the development of monitoring and environmental impact studies, in order to allow mitigating the effects of pollution by detecting early effects on microorganisms and helping to prevent a negative impact on marine biodiversity.

Obtaining the Samples
A total of 28 sediment samples were collected in 2013 and were distributed between 22°58′32.55″S; 70°19′18.72″W; and 23°01′32.97″S; 70°29′50.79″W. Seventeen of these samples were collected by the FDI project with a mini box core of 225 cm 2 at depths of 10 m, 50 m, 70-75 m, and 90-110 m, and 0.5 cm of sediment was extracted from each sample for analysis.These were used to analyze the spatial distribution of abnormal shells within the bay.
The remaining 11 samples belonged to Fondecyt project no.1130511 and were extracted with two sediment cores (ZA and Z1A).These were obtained with polycarbonate gravity core with a size of 68 × 71 mm × 60 cm long, being considered for the analysis of the temporal distribution of foraminifera in the area (Fig. 1).The ZA core was obtained at 75-m depth, and the Z1A at 100-m depth.Fourteen samples belonging to the first 8 cm of each core were analyzed.Subsequently, all the samples were transported to the foraminiferology laboratory of the University of Concepción, where they were dried and sieved into 212-µm, 150µm, and 63-µm fractions.However, only the largest fraction (212 µm) was analyzed, because for this study, we identified only shell deformities, which are defined as abnormalities that are generated throughout the lifetime of each foraminifera [24] and result as a consequence of regeneration following shell damage [25].Smaller fractions (150 µm and 63 µm) were excluded from the analysis, since these malformations refer to abnormalities of ontogenetic development in the juvenile stages of foraminifera, which occur very commonly in early stages of development, even in unpolluted coastal areas [24,26].
To quantify and identify the organisms most efficiently, we subdivided each sample into equal parts with a splitter.The total community was analyzed including living + dead foraminifera.Thus, the succession times differ from each other and could represent specimens that were alive a few years ago or perhaps thousands of years ago [27].Therefore, this could be a limitation of this work.However, we have the dating, and the recent layer is very thin in both cores, so the proportion of living organisms is very low.

Sediment Geochronology
The sediments of the Z1A core were dated (Fig. 2a) using the alpha decay of 210 Po in secular equilibrium with 210 Pb and with a constant rate of supply (CRS) model [29,30], establishing accumulated inventories of 210 Pb xs (not supported) from its exponential decay with depth: where t 210 is the decay constant of 210 Pb (32.18 years), I 0 and I z are the total unsupported inventories of 210 Pb ( 210 Pb xs ) in the sediment core at surface and below depth z, respectively.Inventories (I 0 and I z ) are estimated from: where I is the 210 Pb xs inventory (dpm cm −2 ), Ai is the total 210 Pb xs (dpm g −1 ), ρi is the bulk density in interval i (g cm −3 ), and h is the interval thickness (cm).The sustained activities were deduced from the constant activities measured at the base of the nucleus, in equilibrium with 226 Ra according to exponential decay [29,31].

Heavy Metal Estimation
A surface sample from the dredge was selected.Each sediment sample was previously dried in an oven at a temperature of 105 °C for a period of 2 h.Subsequently,

Identification of Shell Abnormalities
Anomalous shells were separated in the foraminiferology laboratory of the University of Concepción (UdeC) using a stereoscope, and the types of abnormality in each species were classified by sample using an optical microscope.Abnormalities were classified as described by [12,32,33] (Table 1) and were stored in foraminifera holding plates, from which photographic records were taken using a SEM SU3500 Hitachi from the Chillán campus of the UdeC.

Quantitative Analysis and Distribution of Anomalous Shells
The Foraminifera Abnormality Index (FAI) was estimated for each sample.This is defined as the total percentage of abnormal shells [14,34].When the percentages are greater than 1%, they are considered as a community under environmental stress conditions [32,34], being considered as high FAI values.Once the FAI values were obtained, their spatial distribution along the bay (dredge samples) and their temporal distribution through the dated sediments (core   samples) were analyzed.This allowed identifying the locations and time periods that recorded higher proportions of abnormalities.Additionally, the minimum (P min ) and maximum percentages of abnormality (P max ) were calculated for each sample, given that the differentiation between normal and abnormal shells is sometimes confusing.Therefore, P min was estimated as the percentage of strongly abnormal shells, and P max was obtained as the total percentage of abnormality, including those observed with a lower degree of deformity [33,35].P min and P max values less than or equal to 1 were not considered as an effect of environmental stress, as they indicate normal conditions [25,32].

Relationship Between Shell Abnormalities and Heavy Metals
As a first step, the distribution type of the data at the temporal and spatial level was identified with a Shapiro-Wilk normality test, which resulted in all of these data not corresponding to a normal distribution (Annex a).Subsequently, the dependent (types of abnormalities) and independent The GLM model is used as an exploratory analysis from versatile regressions based on the maximum likelihood method [36], in order to evaluate as a hypothesis the positive relationship between the different types of abnormality and the concentrations of heavy metals in sediments.Therefore, a significance value lower than 0.05 (p < 0.05) and the log-likelihood value indicates the best fit of the model.The GLM is usually estimated when: (1) a large number of variables are evaluated.(2) The distribution of the data is not normal.(3) In case of count data (number of shells with a certain type of abnormality) and proportion data (FAIpercentage of abnormality).( 4) When the data has a large number of zeros and, consequently, the variance tends to increase linearly towards the mean [37].
To evaluate the GLM model, a total of 10 samples (core and drag samples) were analyzed, and the N values analyzed for each GLM were between 30 and 278, corresponding to the number of shells with a specific type of abnormality (wrong coiling, reduced chamber size, scars, and number of type of anomalies).These counts were derived from the total number of foraminifera per sample, in which most foraminiferal abundances ranged between 200 and 400 specimens (normal and abnormal individuals), being statistically representative [38,39].Once the significant relationships were obtained from GLM analysis, a principal component analysis (PCA) was carried out with the PAST 4.08 software.It obtained correlations between the abnormalities of wrong coiling, reduced chamber size, scars and number of anomaly types, and the different heavy metals obtained with dredge and sediment core samples.Finally, correlations were made with the load coefficients of the different variables (abnormalities and heavy metals) in the principal components, which explained a greater variability of the PCA, considering a significance of α = 0.05.

Results and Discussion
The results shown below were part of Tavera's master thesis work [40].All the sediment geochronology results were published in Tavera et al. [28].

Sediment Geochronology
The core shows a good exponential decay with very low mixing, which allowed obtaining a good estimate of the sustained activity, in balance with the activity of the parent 226 Ra (Fig. 2b) (2.04 ± 0.71 dpm g −1 ; r 2 = 0.95, p < 0.01).A hiatus at the limit of the excesses (not supported) would seem to correspond to a tsunami event, which, given the ages, would probably correspond to the 1877 event reported for the area.This explains the low activity between Sections 6-7 and 13-14 cm, and being at the limit of the excesses, under this layer it was not possible to estimate the ages with this method.Thus, the excesses are concentrated in the first 6 cm, accounting for a total inventory of 57.8 ± 0.1 dpm cm −2 .The values determined for sustained activity and total inventory would be in the range of what was previously observed for the area (3.4 dpm g −1 -50 dpm cm −2 ) [41].The errors associated with ages increase with depth, which is inherent to the CRS method used.The  estimated sedimentation rate was 0.11 ± 0.02 cm year −1 , similar to other bays with similar conditions in northern Chile [42,43], which allows establishing an approximate stratigraphy of the ZA core.

Heavy Metals in Sediments
The concentrations of metals (Al, Ti, Sr, Ba, Sn, V, Cd, Mn, Cu, Zn, Pb, Mo, and Cr) were estimated from superficial sediment and sediment core samples (Table 2 and Table 3).Most of them could arrive in the marine environment through air pollution (dust) and domestic and industrial sewage but also have important inputs from natural sources; therefore, their accumulation would be related to anthropogenic and natural sources [44].The elemental composition of sediments in Mejillones Bay was in the range of previous measurements in the zone and farther south bays (Table 4).
The Al concentration was more variable than that of Ti, and both elements were dominated by terrigenous sources.Their variability reflects variability in sediment composition and grain size [45].Mn and V are sensitive to redox changes, Mn is precipitated and accumulated in oxic conditions, while V is accumulated in less oxic conditions (reduced sediments), because this element exhibits higher concentrations due to the low oxygen conditions prevailing in Mejillones del Sur Bay.Similarly, Cd showed high concentrations (Table 2 and Table 3) due to its affinity with sulfides in sediments, normally produced during sulfate reduction in anoxic environments [46], and the range of values was higher than that in the southern sites (Inglesa and Caldera bays) because Mejillones sediments are heavily influenced by the OMZ, which helps to maintain reduced conditions.
Upwelling is an important source of Cd [47] and contributes to its accumulation.Cu, Ni, and Zn, elements highly related to the organic flux to the sediments derived from primary productivity [48], showed higher range values, which was more evident for Cu.Most of these elements were influenced by anthropogenic activities in the zone.The concentrations of other elements were not reported previously but were considered in the analysis.Sr is related to the terrigenous supply to the bay sediments as well as waste from industries, such as the disposal of coal ash and incinerator ash.The concentrations reached 3320 mg g −1 in the surface sediments (Table 2).Similarly, Sn and Cr could also be related to industrial and urban waste zones, and their values reached a maximum of 4440 mg g −1 and 157 mg g −1 , respectively.These values could be highly influenced by the industrialization of the zone, because the preindustrial values estimated from core Z1A were 956 mg g −1 and 98 mg g −1 for Sn and Cr, respectively (7-8-cm layer).In the case of Ba, the values were higher in surface sediments compared with pre-industrial values, but its concentration was highly influenced by organic fluxes to the sediments from primary productivity and suboxic conditions [48].

Taxonomic Diversity
All the results related to the taxonomic diversity and abundance of foraminifera are reported in Tavera et al. [28].However, the presence of Bolivinitides, species of Valvulineria and Cassidulina (Annex b), being typical of the Mejillones Bay, was highlighted, which evidences characteristics of hypoxia and with a large amount of organic matter.

Abnormal Shells in Core Samples
The most frequent abnormality type recorded was wrong coiling of the chambers with 39%, followed by reduced chamber sizes (20%) and regeneration marks or scars (18%).On the contrary, shells with aberrant shapes were evident in 7% and shells with double aperture in 5%.The remaining anomaly types (decalcification, double shells, compressed shells, and protuberances in the chamber) all recorded values lower than 4% (Fig. 4a).
Of the 31 species recorded in total, 14 were recorded with abnormal shells.Bolivina seminuda presented the greatest variability and number of abnormalities with 78% Fig. 4 Percentage of a types of abnormalities, b abundance abnormal shells recorded by species, and c types of abnormalities by each layer of ZA and Z1A cores.Percentage of d types of abnormalities, e abundance abnormal shells recorded by species, and f types of abnormalities by each sample obtained with a dredge ◂ distributed in 9 types of abnormality (reduced chamber size, protuberances, wrong coiling, aberrant shape, compressed shells, double shells, additional chambers, scars, and decalcification, Fig. 3), as well as being one of the most predominant species in this area [28].The remaining species showed much lower percentages, with B. plicata making up 4% of the abnormalities, followed by Cassidulina laevigata with 3%.B. striatula and C. variabilis showed even lower percentages (0.56% each) (Fig. 4b).With respect to the total number of abnormalities, wrong coiling was observed in greater proportion in the most superficial layer of the ZA core (46%), while in the 1-2-cm layer of this core, an increase in the proportion of shells with double aperture (19%) and reduced chamber size (26%) was detected (Fig. 4c).

Abnormal Shells in Dredge Samples
Decalcification of the shells was the most frequently registered anomaly (37%), followed by wrong coiling (24%), scars or regeneration marks (12%), and aberrant shapes (10%).On the other hand, shells with reduced chambers, compressed shells, and additional chambers ranged from 4 to 5% (Fig. 4d).These abnormalities were registered in 15 of the 39 species collected in total.Buliminella elegantissima presented the highest proportion of abnormalities (29%), represented by shells with scars and decalcification (Fig. 3).
Consecutively, 19% of the anomalies were evidenced in E. exigua, distributed in six types (reduced chamber size, protuberances, wrong coiling, aberrant shape, scars, and decalcification).Additionally, 16% of anomalies were present in B. seminuda and 8% in Rosalina cora.The remaining species showed between 1 and 5% of abnormalities (Fig. 4e).Regarding the number of total abnormalities per station in the 10-m isobath, there was an increase of 100% of decalcified shells (T4) and 50% with scars (T7).In contrast, wrong coiling increased in proportion at the 50-m isobath stations (Fig. 4f).
This reflects that the bay's hypoxic conditions hinder the mineralization process for the formation of calcium carbonate in the shell.This causes thinner and less ornamented walls [49,50], generating a greater susceptibility to the entry of heavy metals into the cell.Since these elements come from seawater, they are incorporated as benign ions by the cell membrane and are transported intracellularly [50].Once there, they deteriorate the proteins of the cytoskeleton and interfere in calcification, specifically in the formation of calcium crystals.This produces a deformation of the calcium carbonate crystalline structure, generating abnormalities [50,51].
Therefore, the large number of anomalous B. seminuda, B. elegantissima, and E. exigua shells indicates that these species have a morphological response to heavy metal contamination.Therefore, they can be used as contamination bioindicators, a useful tool for environmental impact studies, especially considering their high abundances under hypoxic conditions [28,49] and their easy identification as a species.However, it is important to note that the proportion of anomaly types varied depending on the sampling method, and with it, the species recorded with these morphological effects.This difference may be due to the fact that temporal analyses (core samples) produce a reconstruction of past environmental conditions, where said abnormalities reflect changes in the ecological conditions of the area, such as changes in salinity, organic matter, and oxygen concentration over time [32,33] and therefore these temporal abnormalities do not only indicate the contamination of the bay caused by heavy metals.

Spatial Distribution of the Foraminifera Abnormality Index (FAI)
According to the records obtained by dredge sample, the FAI showed highly variable values spatially but no clear relationship pattern between the FAI index and distance to the coast.In the 10-m isobath, the maximum abnormality value was 6% (station T6), as opposed to that reported in the 50-m isobath which showed a maximum FAI of 16% (station T8).Between 70-and 75-m depth, values of up to 5% were reported (station T7), similar to what was recorded between 90 and 110 m, where the values were 4 and 5% at the 90-m and 100-m stations, respectively (Fig. 5a and Table 5).
The highest FAI values or percentage of abnormality were evidenced off the Angamos Point area (FAI between 4 and 16%); these stations were the furthest away from the thermoelectric power plants and port terminals.The stations close to the Chacaya Point area showed lower proportions of anomalous shells (between 2 and 3%) and even lower in the stations located in the most extreme points of the bay (T1-75 m and BMS-110 m) (Fig. 5a).Therefore, the pattern of currents plays a fundamental role in the distribution of these anomalies, as the interaction of sub-Antarctic and Antarctic water masses and the South Pacific Subtropical Anticyclone (APS) generates a type of cyclonic and anticyclonic eddies inside the bay [52,53].This displaces the water masses in the direction of Angamos Point and accumulates heavy metal concentrations in the sediments to the southwest of the Bay, evidencing high FAI percentages (Fig. 5b and 5e).These values are similar to those described for other polluted areas, registering abnormality rates between 10 and 20% in estuaries in England [15], 3 and 7% in fjords in Norway [32], 2 and 3% off the coasts of Israel [12], and from 3.5 to 19.1% on the coasts of the Adriatic Sea [34].

Temporal Distribution of the Foraminifera Abnormality Index (FAI)
As recorded in the sediment core samples, the highest abnormality indices were recorded in the first 2 cm of the Fig. 6 GLM at temporal distribution samples (ZA and Z1A) with a positive and significant relationship (p < 0.05) between wrong coiling (N = 278 shells) and Ti (a), Mn (b), and Ni (c); reduction in the size of the shell chambers (N = 140 shells) and Ni (d); scars or regen-eration marks (N = 125 shells) and Ti (e) and Mn (f); and number of types of anomalies (N = 30 shells) and V (g).GLM at spatial distribution (drag samples) with a value p < 0.05, obtained between the wrong coiling (N = 88 shells) and Ti (h), Sr (i), and Mn (j) Z1A core, and as the depth of the stratum increased, the FAI was shown to decrease, registering anomalous shells up to 2-3 cm in depth (year 1996).Therefore, from the 4-5cm stratum, no shells with abnormalities were recorded in either of the two sediment samples.Core ZA yielded lower abnormality indices compared to Z1A, ranging from a minimum FAI of 0.2% (1-2 cm) to a maximum of 3% (2-3 cm) (Fig. 5c and Table 5).
Regarding the sediment cores, abnormality indices were recorded in the foraminifera of the most superficial sediments (first 3 cm), corresponding to the years between 1996 and 2012.However, the percentages varied considerably between both cores, reflecting higher proportions in the Z1A core (12%) compared to ZA (3%).Meanwhile, the deeper strata (between 4 and 8 cm) showed an absence of anomalous shells in both sediment cores (Fig. 5c and 5d).
The presence of anomalous shells in the three most superficial centimeters (between the years 1998 and 2012), coincided with the years in which industrial activity was developed in the area (thermal power plants between the years 1996 and 2011 and port terminals in the year 2002).This indicated a morphological response of foraminifera to environmental stress and anthropogenic impact, as described by other authors for contaminated areas [33,[54][55][56].In addition, species that commonly occur in adverse conditions (eutrophication and low oxygen) were dominant, such as Bolivina and Buliminella (Annex b).
In turn, the high FAI values in Z1A compared to ZA may be related to the location of this core in the area near Angamos Point, due to the displacement of sub-Antarctic and Antarctic waters and the APS towards this area [52,53], allowing the accumulation of sediments and heavy metals, as evidenced in other sampling points bordering this location (Fig. 5a).Therefore, the local oceanographic characteristics of the Bay play a fundamental role in the distribution of heavy metals and therefore in abnormal foraminifera shells.This indicated that in the most recent years (between 1996 and 2012), the bay has presented environmental stress conditions generated by anthropogenic impact, which in turn affects the environmental conditions of the area (eutrophication and microxic conditions).This has been described for other studies in Mejillones Bay, where negative effects on the benthic macroinvertebrate community [22] and bioaccumulation of heavy metals from the crustacean Emerita analoga [23] have been found.

Minimum (P min ) and Maximum (P max ) Percentages of Abnormality
In the samples obtained with the dredge, the highest percentages were evident at a depth of 50 m, where the P min was 6% and the P max was 12%, followed by the 90-m isobath which presented values of 9% for both P min and P max .On the other hand, in the foraminifera collected with sediment cores, the 0-1-cm stratum of the Z1A core recorded the highest values (P max of 11% and P min of 4%), and consecutively the 1-2-cm stratum of this same core showed 4% and 3% of P max and P min , respectively (Fig. 5e).
The P min and P max values found were within the ranges described under environmental stress conditions, similar to what was reported in estuaries and polluted lagoons, registering P min values between 2 and 21% and a P max values between 2 and 27% [33,35].However, some values less than or equal to 1% were recorded (Fig. 5e), which were therefore not considered environmental stress effects as described by Alve [32] and Stouff et al. [25].

Relationship Between Shell Abnormalities and Heavy Metals
According to results obtained by the generalized linear model (GLM) at temporal distribution samples (sediment cores), positive and significant relationships (p < 0.05) were identified between wrong coiling of the chambers and Ti, Mn, and Ni concentrations (Fig. 6a-c).Also, a significant relationship between reduced chamber sizes and Ni levels (Fig. 6d), and the scars in the shells with Ti and Mn concentrations were established (Fig. 6e-f).Additionally, this same positive and significant relationship was estimated between the number of types of morphological anomalies with vanadium (V) (Fig. 6g).On the other hand, in the spatial distribution samples (dredge), the GLM yielded a positive and significative relationship (p < 0.05) between wrong coiling with Ti, Sr, and Mn levels (Fig. 6h-j).Therefore, the hypothesis proposed for this model at the temporal and spatial level is accepted, finding a significant relationship between several types of abnormalities (wrong coiling, reduced chamber size, and scars) with some heavy metals (Ti, Mn, Ni, V, and Sr).
These GLM results were corroborated by principal component analysis (PCA), performed for the sediment core samples.This indicated a positive correlation between the previously described response variables (wrong coiling, decreased chamber size) with Ti, V, Mn, and Ni concentrations (principal component 1).On the other hand, foraminifera with scars were positively and significantly correlated only with Sr concentration (principal component 2).We believe that an increase of Sr in the sediment could displace Ca in the shells and be the cause of this deformation (Fig. 7a, 7b and Annex c and d).In contrast, the results obtained from the PCA of the spatially distributed samples (drag samples) did not agree with the spatial GLM (Fig. 6h-j), because the correlations that showed positive significance corresponded to the number of anomaly types with Ti and Ba concentrations (principal component 2) (Fig. 7b and 7c, Annex c and d).
Despite the fact that several previous studies have described a relationship between the increase of Cd, Zn, Pb, Cu, and Hg and abnormal foraminifera shells (e.g., [57][58][59]), this was not found in our analyses, in which a positive and significant correlation was evidenced between Ti, V, Mn, Ni, and Ba and different types of abnormalities (wrong coiling, scars, decreased chamber size) and number of abnormality types, but not specifically with the FAI.This is the first work that describes these relationships.This discrepancy in the identification of the metals that generate shell abnormalities may be related to their bioavailability in the sediments, which varies depending on the environmental, geochemical, and biological factors of the area [56,60], as well as the enrichment of the sediments in said metals, since the development of morphological anomalies is evidenced when the concentration of heavy metals accumulated in the shell exceeds the thresholds supported for each species [61].Therefore, in this case, the effect of Mn on abnormalities specifically reflects an enrichment of Mn ions in the most recent sediments (0-1-cm core Z1A, Table 1), which generates the replacement of commonly used Mg ions in the generation of porcelaneous foraminifera shells, indicating toxicity effects [55,62], as reflected in the abnormalities of species of the genus Quinqueloculina (Fig. 3l).
Similarly, the enrichment of Ni in the most recent sediments of the bay (123 and 171 mg g −1 , Table 2) favored the positive and significant correlations with various types of anomalies (wrong coiling, decreased chamber size, and scars).The same was described for other contaminated areas on the island of Mayorca, where Ni was found in the sediment between 0.83 and 4.14 mg g −1 , correlating only with protuberances and abnormal chamber growth [63].This metal comes mostly from anthropogenic activities and enters the coast through terrigenous material contribution [64].On the other hand, despite the fact that V has been incorporated in various studies (e.g.[16,51,65],), it had not been found to be correlated with shell abnormalities.However, the positive and significant correlation with various types of anomalies in this study was given by the affinity of V for redox environments generated in anoxic zones, typical conditions in the Mejillones Bay, where sinkholes and sediment enrichment with this metal occur, facilitating the incorporation of these ions into foraminifera shells [66], affecting the calcification process.On the other hand, deformities have been found in foraminifera caused by Ti concentrations in contaminated sediments (up to 465 mg g −1 ), which tends to accumulate in silt-clay substrates [50].This is similar to the Mejillones Bay, evidencing this same type of sediments and with maximum concentrations of 2570 mg g −1 (Table 1).

Conclusions
The species B. seminuda, B. elegantissima, and E. exigua presented a greater susceptibility to environmental stress in Mejillones Bay, evidencing a greater number of shell anomalies, considered early biomarkers of heavy metal contamination.On the other hand, variation in the percentages and types of abnormalities in foraminifera in the Bay indicated changes in environmental conditions spatially and over time.
Although these are the first results of abnormalities in foraminifera obtained in Mejillones Bay and the first study carried out in the coasts of Chile, these results coincide with those found in other polluted coastal areas (e.g., estuaries of England, coasts of Israel, and fjords of Norway).In addition, negative effects on community diversity and bioaccumulation of heavy metals in benthic invertebrates have been previously reported in this bay.
The spatial distribution of the FAI and of the heavy metals in the sediments depends on the oceanographic conditions of the area, where the cyclonic and anticyclonic eddies inside the bay move the water masses in the direction of Angamos Point so abnormalities and metals accumulate southwest of the Bay, while, at a temporal level, the FAI values were evident in the most recent years (between 1996 and 2012) and were found in an even greater proportion in the last year (2012), coinciding with the period in which industrial development began in the bay.This indicated environmental stress conditions in the area.
Despite not finding a correlation between the FAI and heavy metals in the sediments, significant specific correlations were identified between Ti, Mn, Ni, V, and Ba with wrong coiling, scars, decreased chamber size, and number of anomaly types, different from what was found in previous studies.This evidenced that there is no generality of heavy metals that cause toxic effects in foraminifera, and therefore, their adverse effects depend on the particular geochemical conditions of the area (sediment enrichment, redox conditions, affinity with the type of substrate), which contribute to a greater bioavailability of heavy metals and replacement of ions necessary for shell generation.

( 1 )
tz = t 210 * ln I z ∕I 0 (2) I = ∑ Ai ih the samples were ground in an agate mortar and weighed on an analytical scale.Heavy metal concentrations were obtained by wavelength X-ray fluorescence (WDXRF) using an ARL ADVANT`XP + sequential model from THERMO, where the samples were analyzed in a helium environment.

Fig. 3 a
Fig. 3 a B. seminuda shell without abnormalities.Types of abnormalities recorded in B. seminuda shells: b scars, c reduced chamber sizes, d decalcification, e double aperture, f, g additional chambers.h B. elegantissima shell without abnormalities, i scars, and j decalcification in B. elegantissima.k Enlarged opening in H. depressula, l scars in Q. seminula, m scars in T. gramen, and n additional chambers in N.

Fig. 5 a
Fig. 5 a Spatial distribution of the Foraminifera Abnormality Index (FAI).b FAI values (%) estimated for samples obtained with a dredge.c Temporal and vertical distribution of the Foraminifera Abnormality Index (FAI).d FAI values (%) estimated for the samples

Fig. 7
Fig. 7 Principal component analysis (PCA) and variance (%) for each component in temporal distribution samples, considering morphological anomalies (a) and heavy metals (b).PCA and variance (%) for

Table 1 Types of abnormalities in benthic foraminifera shells reported by authors and considered for this study
[28]ndicates industrial activity in Mejillones Bay.Extracted from Tavera et al.[28]

Table 2
Trace metals (mg g −1 ) measured in core ZA and Z1A at Mejillones Bay.sd, standard deviation from three replicates; nm, not measured variables (Al, Ti, Sr, Ba, I, Sn, Cl, V, Cd, Mn, Cu, Zr, Zn, Pb, Mo, Cr, Ni, Rb, W, Cs, Co, As, Ta, Nb, Pt, Os, and Y) were identified.Once the variables were identified, a separate generalized linear model (GLM) was performed for the samples at the temporal (layers cores) and spatial (dredge samples) level, which presented abnormal shells.

Table 4
Range values of metal concentrations (min-max) in sediments under similar oceanographic conditions in northern Chilean continental margin

Table 5
Number of abnormal shells, total foraminifera abundance, and FAI values from dredge samples (spatial distribution)