Sabrina Tartu
email: tartu.sabrina@gmail.com

Pierre Blévin

Paco Bustamante

Frédéric Angelier

Claus Bech

Jan Ove Bustnes

Melissa Chierici

Agneta Fransson

Geir W Gabrielsen

Aurélie Goutte

Børge Moe

Christophe Sauser

Julien Sire

Christophe Barbraud

Olivier Chastel

A U-Turn for Mercury Concentrations over 20 Years: How Do Environmental Conditions Affect Exposure in Arctic Seabirds?

Keywords: Global changes, Trace elements, Temporal trends, Rissa tridactyla, Chlorophyll a, Sea-ice extent, Feeding ecology

published or not. The documents may come

TOC Art

Synopsis: Prey community shifts and MeHg release from young sea-ice are counterbalancing the previously observed decline of Hg in kittiwakes.

INTRODUCTION

The Arctic is warming two to three times faster than any other region on Earth 1,2 . Key indicators such as temperature, precipitation, snow cover, sea-ice thickness and extent, show rapid and widespread changes underway in the Arctic. Between 1971 and 2019, the annual surface temperature changes in the Arctic (land and ocean) were three times higher than the global average during the same period 3 . The loss of sea-ice cover over extensive areas of the Arctic Ocean has allowed radiation to penetrate the water during summer, storing heat in the upper ocean. Unique ecosystems, such as those associated with multi-year sea-ice or millennia-old ice shelves, are at risk and some are vanishing 3 . These rapid changes are causing fundamental changes in ecosystems that affect biogeochemical fluxes, bottom-up processes, ecosystems, and food-webs, which may lead to modifications in contaminant exposure in the Arctic biota 5,6 .

Among contaminants of high concern, mercury (Hg) is a non-essential trace element from both natural and anthropogenic sources. The Arctic has been considered as a sink for Hg deposition, with the highest deposition observed in the European part of the Arctic, downwind of major anthropogenic source regions 7 . However, the fate of Hg trends in Arctic ecosystems in light of global change remains unclear.

Global warming, which will lead to an earlier onset of thawing and a later start of freezing, is likely to extend the season of Hg methylation 4 . Methylmercury (MeHg), which is known to be the most toxic among the Hg compounds, is generated when inorganic Hg is dissolved into water 8 . While bacteria methylation is considered as the major pathway 9 , chemical reactions yielded by environmental factors such as pH and temperature contribute to make MeHg more bioavailable in aquatic ecosystems 10 . MeHg has been identified as a priority contaminant to consider in future climate change-pollutant scenarios 11,12 . An increased period during which MeHg can be produced may then lead to enhanced MeHg exposure in marine food webs. In contrast, global warming is also likely to modify Hg demethylation processes in polar ecosystems. One of the most important degradation processes of MeHg in ocean surfaces is photo-induced demethylation. MeHg decomposes under ultraviolet and visible light 14 .

Therefore, a reduced sea-ice extent will expose a larger ocean surface and could enhance Hg photo-demethylation. The transfers and transformations of Hg in Arctic aquatic systems are poorly understood and need more research 15 . Hg concentrations biomagnify throughout the food webs, with higher concentrations in higher trophic level species 16 . Modifications in productivity at the lower levels of food webs may affect the amount of energy available for predators to grow and reproduce. Energy intake, as measured by body condition or body size, and feeding habits are therefore interesting variables to consider when studying Hg concentrations in higher trophic level species.

The dynamics of Hg exposure over time can be monitored through the exposure of top predators that are considered as relevant bioindicators of environmental pollution [17][18][19] . A review paper has suggested that, as a consequence of global warming, the temporal trends of MeHg in top predators are uncertain 12 . Increasing trends of Hg concentrations were observed in several Arctic seabird species 20,21 , and have been attributed to anthropogenic emissions, climate-related variables, or changes in food web contamination. The Barents Sea is a hotspot of Arctic warming. Between 1979 and 2013, the number of sea-ice free days around the Svalbard Archipelago (Norwegian Arctic) has increased by 41.8 days, which is far greater than what had been monitored in other parts of the Arctic during the same period 22 . Since the 1990s, Hg exposure followed a non-linear increase in some mammals of Svalbard and have been attributed to changes in the food web, as well as re-emission of previously stored Hg from thawing seaice, glaciers, and permafrost 18,19 .

Black-legged kittiwake Rissa tridactyla (hereafter 'kittiwakes'), are generalist top predators, and can be useful indicators of Hg contamination in their habitat 23 . Between 1975 and 2010, at the exception of a few 24 , kittiwake colonies have been declining throughout their breeding range, and this coincided with an abrupt warming of sea-surface temperature in the 1990s 25. In addition, in this species, high Hg concentrations have been associated to decreased breeding probability, lower concentrations of the parental hormone prolactin, and decreased breeding success 23,26,27 . Although research conducted on kittiwakes was based on field studies, several experimental studies conducted on other bird species have led to similar conclusions 28,29 . For instance, Hg exposure led to lower survival and hatching success in bird laboratory models 28,29 . Therefore, the observed decline in kittiwake colony size 25 may be aggravated by Hg-related breeding failure 30 .

Kongsfjorden is a fjord located on the west coast of Svalbard, where kittiwakes breed in large colonies. In recent years, Kongsfjorden has been strongly influenced by the West Spitsbergen Current, with increasing temperatures and declining sea-ice cover inshore. Since 2005, an increasing proportion of Atlantic species has been part of kittiwakes diet 31 , which may possibly influence Hg exposure. The aim of the present study was to investigate the temporal trends of blood Hg exposure in kittiwakes over two decades and identify its potential drivers. We predicted that, in concordance to what is observed in other seabird species from Canada and Greenland and in Svalbard mammals 18,19,21 , a climate-driven diet modification would lead to a non-linear increase of Hg in Svalbard kittiwakes. In addition, since fish community modifications in Kongsforden were driven by sea-ice extent 31 , we predict that sea-ice extent, over other climate-related variables, would be the strongest predictor of Hg concentrations in Svalbard kittiwakes.

MATERIAL AND METHODS

Fieldwork and blood analyses

Captures, blood sampling, and analyses are described in details in the supplementary information. Briefly, from 2000 to 2019 (excluding 2010), 1152 kittiwakes were captured during different phases of the breeding cycle in Kongsfjorden, Svalbard (78°54'N, 12°13'E), Norway. A blood sample (ca. 1 mL) was collected from the alar vein. Birds were weighed to the nearest 5 g using a Pesola spring balance, and their skull length was measured to the nearest 0.5 mm with a sliding calliper. Birds were captured during the breeding period, which can be separated into three stages: the pre-laying, the egg-incubation, and the chick-rearing phases.

Temporal trends of total blood Hg (hence 'Hg') and associations with climate-related variables were analysed in chick-rearing kittiwakes only (n = 392), in order to reflect local Hg exposure and to make between-year comparisons possible. Methods for QA/QC, molecular sexing, Hg, and stable isotope analyses are detailed in the supplementary information. An overview of the samples available per year and per breeding stage is given in Table S1.

Environmental and biological data

Mercury, carbon and nitrogen stable isotopes in blood reflect the diet during the previous month, up to the two previous months [32][33][34] , see supplementary information for more details. Thus, sampling conducted during the chick-rearing stage represents the feeding ecology within the two previous months (May-June). We selected 10 environmental variables known to influence Hg in the ecological niche of kittiwakes 4,21 . These variables included the North Atlantic Oscillation (NAO), the Arctic Oscillation (AO), atmospheric gaseous elemental Hg (GEM), air temperature, sea-ice extent, sea surface temperature anomaly (SSTA), chlorophyll a concentrations (Chl a) as a proxy of primary productivity and a useful indicator of fisheries production at both regional and global scales 35,36 , and precipitation. We also analysed the relationships between Hg concentrations and ocean acidification indicators such as pH and aragonite saturation at the bottom of the water column (see supplementary information).

Except for global measures (NAO and AO), all measures were taken locally during the spring period in the foraging area of kittiwakes. Climate-related variables, sampling periods, and sampling years available for each variable are summarized in Table S2. Except for pH and aragonite saturation which will be the subject of another paper, temporal trends of the eight other variables suspected to predict Hg concentrations in kittiwakes are shown in Figure S1.

Statistical analyses

We detrended biological-and climate-related time-series in order to test if year-to-year variations in the response variable are driven by year-to-year variations in some explanatory variables. To do so we extracted the residuals from independent generalized additive mixedeffect models (GAMMs) with year as the predictor, with a smoothed parameter k adapted to the shape of the trend (e.g. k=4 for Hg and climate-related variables, k=5 for morphometrics and k=8 for stable isotopes). Then we examined the residual covariation between Hg and biotic or abiotic variables 37 . All analyses were conducted using detrended variables except those aimed specifically at testing the effect of year on a variable (e.g. temporal trends of Hg). We provide all statistical methods in details in addition to results conducted on non-detrended data in the supplementary information.

We first tested Hg and stable isotopes variations (detrended variables) across the breeding cycle by using linear mixed-effect models (LMMs) with breeding period and sex as independent categorical 3-and 2-level predictors, and bird identity as a random factor. Then, we investigated whether sex, body mass, size corrected mass SCM (i.e. residuals from the mass ~ skull size regression), scaled mass index SMI 38 and skull size would influence individual detrended Hg concentrations using LMMs. We obtained estimated marginal means (EMMs) for specified factor combinations of the LMM (breeding period × sex interactions) by using the 'emmeans' package 39 .

Second, to investigate the temporal trends of Hg and other biological variables (δ 13 C, δ 15 N, body mass and skull size) in chick-rearing kittiwakes, over a 20-year period, we used GAMMs on non-detrended data to account for non-linear relationships between response variables and predictors, and for recaptured individuals. Year, δ 13 C or δ 15 N and sex were used as continuous predictors, and bird identity as a random factor to analyse temporal trends of Hg, whereas models explaining the temporal trends of δ 13 C, δ 15 N, body mass and skull size only used sex as a predictor and bird identity as a random factor. For each model, we adapted the smoothed term (k) to avoid overfitting. Break-points were identified by using a Bayesian approach ('broken stick' modelling technique 40,[START_REF] Thiebot | Gentlemen First? 'Broken Stick' Modelling Reveals Sex-Related Homing Decision Date in Migrating Seabirds[END_REF] detailed in supplementary information.

Third, we analysed the relationships between detrended Hg concentrations and detrended climate-related variables. To avoid pseudoreplication due to the large, and sometimes unbalanced, number of kittiwakes caught per year, we calculated for each year the average Hg concentrations, δ 13 C and δ 15 N values for each sex, and performed statistical analyses on this subset of n=38 sampling units. Since Hg trends in biota often present non-linear variations [START_REF] Gaden | Mercury Trends in Ringed Seals (Phoca Hispida) from the Western Canadian Arctic since 1973: Associations with Length of Ice-Free Season[END_REF][START_REF] Rigét | Temporal Trends of Hg in Arctic Biota, an Update[END_REF][START_REF] Zabala | Early Breeding Failure in Birds Due to Environmental Toxins: A Potentially Powerful but Hidden Effect of Contamination[END_REF] and have been increasing non-linearly in Svalbard vertebrates 18,19 , we tested both linear (e.g. sea-ice extent) and quadratic (e.g. sea-ice extent + sea-ice extent 2 ) relationships using linear models (LMs). We computed a list of 34 biologically relevant model candidates by adding the variables selected from the best model from the analysis described on step two (Hg variations in relation to morphometrics and sex), to which we added one of the eight climate-related variables without missing years, either as linear or quadratic. We then selected the best fits (linear or quadratic) for all variables included in the top 10 models to explain detrended Hg concentrations, and used conditional model averaging to make inference from all of the models [START_REF] Burnham | Model Selection and Multimodel Inference: A Practical Information-Theoretic Approach[END_REF][START_REF] Lukacs | Model Selection Bias and Freedman's Paradox[END_REF] . Variables with missing years (pH and aragonite saturation at the bottom of the water column) were analysed separately using LMs. Finally, we tested which of the strongest climate-related variables would better explain feeding habits (δ 13 C or δ 15 N), to do so we used LMs with detrended variables.

All models were ranked using AIC [START_REF] Akaike | Factor Analysis and AIC[END_REF] . We obtained the degrees of freedom (df), the log likelihood (logLik), the difference in AIC values between the "best" model and the model at hand (ΔAIC), and a normalized weight of evidence in favor of the specific model, relative to the whole set of candidate models, derived by e (-0.5(ΔAICc)) (weights). We obtained estimates and 95% confidence intervals (95% CI) for each parameter, and an effect was considered significant when the 95% CI did not include 0. Model fit was assessed by using residual diagnostic plots.

All analyses and figures were performed using R software v.3.6.3 [START_REF] Core | R: A Language and Environment for Statistical Computing[END_REF] with the following packages: car, mgcv, nlme, effects, ggplot2 and corrplot [START_REF] Fox | An {R} Companion to Applied Regression[END_REF][START_REF] Fox | Visualizing Fit and Lack of Fit in Complex Regression Models with Predictor Effect Plots and Partial Residuals[END_REF][START_REF] Pinheiro | Nlme: Linear and Nonlinear Mixed Effects Models[END_REF][START_REF] Wei | Corrplot": Visualization of a Correlation Matrix[END_REF][53][START_REF] Wood | Fast Stable Restricted Maximum Likelihood and Marginal Likelihood Estimation of Semiparametric Generalized Linear Models[END_REF] . Data sharing policy is detailed in supplementary information.

RESULTS

Biological variables influencing individual Hg concentrations

From 2000 to 2019, Hg was analysed in 1152 kittiwake blood samples. Among these birds, biometrics (skull and body mass) were available for 1100 birds and stable isotopes were analysed in 508 birds (Table S1). Mercury concentrations decreased by 55.7% across the breeding cycle, ranging from 1.97±0.03 μg/g dw in pre-laying kittiwakes to 1.09±0.02 μg/g dw (values represent raw mean ± standard error concentrations) in chick-rearing kittiwakes (Figure 1A, Table 1A). This decrease of Hg across the breeding cycle was further supported by analyses conducted on a subset of birds, which were recaptured during different breeding phases within a year (Figure S2, Table 1A). Blood Hg concentrations in males were higher by 0.22 µg/g (0.04; 0.39) during the pre-laying period, 0.47 µg/g (0.33; 0.66) during the incubation period and 0.37 µg/g (0.21; 0.44) during the chick-rearing period, in comparison to females (values are estimates from EMMs and 95% CIs in parentheses). In 2011, which was the only year where stable isotopes were analysed in a large enough number of individuals for each period, we analysed δ 13 C and δ 15 N variations across the breeding cycle (Figure 1B). Stable isotope values in incubating and chick-rearing birds were lower in comparison to values in pre-laying birds (Table 1B, Figure 1B). In 2011, stable isotopes were not related to sex at any stage of the breeding season (males vs females, δ 15 N: pre-laying period: -0.09 (-0.54; 0.36), incubation: 0.32 (-0.21; 0.85), chick-rearing 0.29 (-0.29; 0.86); δ 13 C: pre-laying period: 0.08 (-0.14; 0.31), incubation: -0.0003 (-0.25; 0.25), chick-rearing -0.15 (-0.42; 0.12), EMMs estimates and 95% CI).

When analysing the effects of sex, body metrics and condition on kittiwake's Hg concentrations, sex over morphometrics explained Hg concentrations (Table 1), we therefore included sex as a fixed factor in further models. Among body metrics and condition calculations, body mass, skull size and SCM were positively related to detrended Hg concentrations (Table 1). However, those relationships were sex-and breeding perioddependent. For instance, in females Hg and body mass were unrelated in pre-laying birds, negatively related in incubating birds and positively related in chick-rearing birds; whereas for males Hg and body mass were unrelated when considering breeding periods separately (Table S3). For all Hg-biometrics relationships according to sex and breeding period see Table S3.

Temporal trends of blood Hg concentrations in chick-rearing kittiwakes

In chick-rearing kittiwakes (n=392), the best ranked model to explain temporal variations of blood Hg included year, sex, and δ 15 N (Figure 2, Table 2). Blood Hg concentrations followed a U-shaped relationship according to year, with an identified change-point in 2014 (Figure 2).

Mercury concentrations decreased by 3% per year until 2013 (LMM: -0.03 (-0.04; -0.02)) and increased by 11% per year after 2014 (LMM: 0.11 (0.06; 0.15)).

Over the 20-year period, δ 15 N and δ 13 C presented oscillating variations over time (Figure 3A,C). Overall, δ 13 C values presented a decrease over time (LMM: -0.019 (-0.023; -0.010)), mostly related to δ 13 C values which remained low after 2015 (Figure 3A), whereas δ 15 N values presented a drop from 2013 to 2015, but no general trend overtime (LMM: -0.004 (-0.017; 0.009), Figure 3C). Body mass and skull size temporal trends were bimodal, with the highest values in 2004 and 2014 for body mass (Figure 3B), and 2004 and 2016 for skull length (Figure 3D). While body mass did not follow a general trend overtime (LMM: 0.15 (-0.38; 0.68), Figure 3E), skull length overall increased by 7% per year (LMM: 0.07 (0.03; 0.10)) (Figure 3D). For Hg, δ 15 N, δ 13 C, body mass, and skull size, we give estimates and 95% CI for the fixed terms and effective degrees of freedom of the smoothers in Table S4 and diagnostic plots in Figures S3,S4.

Environmental variables explaining Hg concentration in kittiwakes

Among climate-related variables, we show the ten most competitive models in Table S5 for both detrended and raw data. The best predictor of detrended Hg concentrations was Chl a (Table 3). Sea-ice extent was the top model when considering non-detrended variables (Table S5, S6). Hg concentrations followed a U-shaped relationship with Chl a (Figure 4) and sea-ice extent. Nitrogen stable isotopes in blood were better explained by Chl a, sex and precipitation, whereas carbon stable isotopes were better explained by Chl a and precipitation (Table 4).

Nitrogen and carbon stable isotope ratios were negatively related to Chl a (Figure 4, Table 4).

Ocean acidification measures, such as pH and aragonite saturation at the bottom of the water column in Kongsfjorden were unrelated to detrended Hg concentrations in kittiwakes

DISCUSSION

The present study documents the temporal trends of blood Hg in chick-rearing kittiwakes sampled in Svalbard (2000-2019). Hg showed a U-shaped temporal trend over the 20-year period, declining by 3% per year from 2000-2013, then increasing by 11% per year from 2014 to 2019At the organismal level, we show that sex is a better predictor of Hg concentrations in comparison to body size or body condition. Hg was higher in males than in females and declined through the breeding season. At the ecosystem level, Chl a concentration predicted Hg concentrations, probably through an effect on kittiwake diet as depicted by the strong relationships between Chl a and feeding habits gauged by stable isotope values.

Biological variables influencing individual Hg concentrations

Hg concentrations in kittiwakes were better predicted by sex. In birds, three ancillary hypotheses are proposed to explain higher Hg concentrations in males than in females: maternal transfer to eggs, body-size dimorphism (i.e. the larger sex consumes larger prey which are more contaminated) and foraging guild differences [START_REF] Robinson | Sex Differences in Mercury Contamination of Birds: Testing Multiple Hypotheses with Meta-Analysis[END_REF] . In kittiwakes, those three hypotheses could partly explain the higher Hg concentrations observed in males, as they are larger and at least during the chick-rearing period males feed at a higher trophic level and on more oceanic prey, as suggested by the higher nitrogen and carbon stable isotope values. However, the higher Hg concentrations observed in males could also result from sex-related physiological differences in Hg retention/excretion [START_REF] Burger | A Framework and Methods for Incorporating Gender-Related Issues in Wildlife Risk Assessment: Gender-Related Differences in Metal Levels and Other Contaminants as a Case Study[END_REF][START_REF] Burger | Methodologies, Bioindicators, and Biomarkers for Assessing Gender-Related Differences in Wildlife Exposed to Environmental Chemicals[END_REF][START_REF] Vahter | Gender Differences in the Disposition and Toxicity of Metals[END_REF] . An experimental study has reported higher fat mass in birds exposed to higher concentrations of Hg [START_REF] Gerson | Environmentally Relevant Methylmercury Exposure Reduces the Metabolic Scope of a Model Songbird[END_REF] . Moreover, in rodents a greater urinary excretion of Hg has been observed in females [START_REF] Hazelhoff | Renal Expression of Organic Anion Transporters Is Modified after Mercuric Chloride Exposure: Gender-Related Differences[END_REF] and despite similar Hg intake, males attained a 20-30% significantly higher whole-body retention in males than females [START_REF] Ekstrand | Mercury Toxicokinetics-Dependency on Strain and Gender[END_REF] .

As previously observed in kittiwakes sampled in Kongsfjorden [START_REF] Øverjordet | Toxic and Essential Elements Changed in Black-Legged Kittiwakes (Rissa Tridactyla) during Their Stay in an Arctic Breeding Area[END_REF] , Hg concentrations decreased as the breeding cycle progressed, and this was not due to sampling bias, as evidenced by repeated captures of the same individual within a breeding season (Figure S4). This decline was observed in females and males, which suggests that female excretion of Hg into eggs does not fully explain the temporal within-season trend. The continuous decrease of Hg is probably related to temporal shifts in prey over the breeding season, coupled to blood Hg sequestration in protein-rich tissues like muscles [START_REF] Peterson | Foraging and Fasting Can Influence Contaminant Concentrations in Animals: An Example with Mercury Contamination in a Free-Ranging Marine Mammal[END_REF] . Since Hg concentrations in red blood cells reflect the exposure of the Hg ingested for a period of up to two months 33,34 , blood Hg concentrations measured during the pre-laying phase partly correspond to exposure in the wintering ground.

Kittiwakes probably feed on higher trophic level more Hg contaminated prey during winter, as suggested by δ 15 N values, although we cannot exclude a Hg and/or δ 15 N baseline shift between the wintering and the breeding grounds. Still, it should be noted that kittiwakes are exposed to higher levels of Hg on their wintering grounds in comparison to that in their reproduction area, as observed in other Arctic seabird species [START_REF] Albert | Seasonal Variation of Mercury Contamination in Arctic Seabirds: A Pan-Arctic Assessment[END_REF][START_REF] Fleishman | Wintering in the Western Subarctic Pacific Increases Mercury Contamination of Red-Legged Kittiwakes[END_REF] ; and those levels are high enough to disrupt reproduction and parental care in this species 23,[START_REF] Tartu | Mercury Exposure, Stress and Prolactin Secretion in an Arctic Seabird: An Experimental Study[END_REF] . This raises a concern about Hg contamination in North Atlantic waters, the area west of the Mid-Atlantic Ridge, has been estimated as a wintering area to 80% of the 4.5 million Atlantic adult kittiwakes [START_REF] Frederiksen | Multicolony Tracking Reveals the Winter Distribution of a Pelagic Seabird on an Ocean Basin Scale[END_REF] and has already been characterised as a hotspot for Hg contamination [START_REF] Fort | Spatial Ecotoxicology: Migratory Arctic Seabirds Are Exposed to Mercury Contamination While Overwintering in the Northwest Atlantic[END_REF][START_REF] Renedo | Contrasting Spatial and Seasonal Trends of Methylmercury Exposure Pathways of Arctic Seabirds: Combination of Large-Scale Tracking and Stable Isotopic Approaches[END_REF] with a likely main source of Hg provided from the Arctic ocean [START_REF] Petrova | Mercury Species Export from the Arctic to the Atlantic Ocean[END_REF] .

Temporal trends of Hg in Kongsfjorden kittiwakes

In adult kittiwakes sampled during the chick-rearing period in Kongsfjorden (Svalbard), blood Hg concentrations decreased by 3 % per year from 2000 to 2013, and increased by 11% per year from 2014 to 2019. Stable isotope values showed oscillating trends over the study period.

We observed a general drop in the values of both δ 15 N and δ 13 C values around 2013-2015, which could coincide with the rise in Hg concentrations in the blood. Since the 2000s, there has been an expansion of the distributional ranges of many northern Atlantic fish communities into the Arctic 31,[START_REF] Fossheim | Recent Warming Leads to a Rapid Borealization of Fish Communities in the Arctic[END_REF] . The increasing presence of these communities coincides with an increase in the intrusion of Atlantic water with associated Atlantic prey into Svalbard [START_REF] Fossheim | Recent Warming Leads to a Rapid Borealization of Fish Communities in the Arctic[END_REF][START_REF] Pavlov | Warming of Atlantic Water in Two West Spitsbergen Fjords over the Last Century (1912-2009)[END_REF] . Periodic arrivals of specific species could result in the yearly oscillations of δ 15 N and δ 13 C values in kittiwake blood. The reason for the oscillating pattern of nitrogen and carbon stable isotope temporal trends, could result from the fact that stable isotope baseline levels can vary globally and might have more diverse sources than Hg [START_REF] Sackett | Nitrogen, and Mercury Isotope Evidence for the Biogeochemical History of Mercury in Hawaiian Marine Bottomfish[END_REF][START_REF] De La Vega | Temporal and Spatial Trends in Marine Carbon Isotopes in the Arctic Ocean and Implications for Food Web Studies[END_REF][START_REF] Pethybridge | A Global Meta-Analysis of Marine Predator Nitrogen Stable Isotopes: Relationships between Trophic Structure and Environmental Conditions[END_REF] . Studying the stable isotopes of Hg could help better track the origin of Hg contamination in kittiwakes [START_REF] Li | Environmental Origins of Methylmercury Accumulated in Subarctic Estuarine Fish Indicated by Mercury Stable Isotopes[END_REF][START_REF] Renedo | Identification of Sources and Bioaccumulation Pathways of MeHg in Subantarctic Penguins: A Stable Isotopic Investigation[END_REF] . The inclusion of sulphur stable isotopes (δ 34 S) in addition to δ 13 C and δ 15 N in further analysis, could provide a better understanding of Hg variations, as δ 34 S has been described as a reliable trophic marker to predict Hg concentrations in six seabird species from the Pacific Canada [START_REF] Elliott | Origin of Sulfur in Diet Drives Spatial and Temporal Mercury Trends in Seabird Eggs From Pacific Canada 1968-2015[END_REF] .

Whereas body mass does not follow a particular trend, body size shows a slight increase over time. This relationship could also result from a change in prey community as suggested by a previous study [START_REF] Christensen-Dalsgaard | Prevailing Weather Conditions and Diet Composition Affect Chick Growth and Survival in the Black-Legged Kittiwake[END_REF] . For instance, kittiwake chicks from a Norwegian colony grew faster in years with a higher proportion of sandeel Ammodytes spp. in their diet [START_REF] Christensen-Dalsgaard | Prevailing Weather Conditions and Diet Composition Affect Chick Growth and Survival in the Black-Legged Kittiwake[END_REF] , they were therefore larger at fledging. Temporal trends of Hg in Arctic seabird communities are inconsistent. Mercury time-series obtained in Canadian Arctic seabirds (eggs and liver) from the 1970s to the 2010s show an increase, followed by a plateauing around 2000 then a decline [START_REF] Braune | Temporal Trends of Mercury in Eggs of Five Sympatrically Breeding Seabird Species in the Canadian Arctic[END_REF][START_REF] Mallory | Do Concentrations in Eggs and Liver Tissue Tell the Same Story of Temporal Trends of Mercury in High Arctic Seabirds?[END_REF] . In contrast, an increase in feather Hg concentrations has been observed for little auks from East Greenland 82 , ivory gulls Pagophila eburnea from the Canadian Arctic and western Greenland 20 , and in eggs of northern fulmars Fulmarus glacialis and Brunnich's guillemot Uria lomvia from the Canadian High Arctic 21 , for the periods of 2007-2014, 1877-2007 and 1975-2014, respectively. These discrepancies likely result from Hg integration time in the sampled tissue [START_REF] Braune | Variation in Organochlorine and Mercury Levels in First and Replacement Eggs of a Single-Egg Clutch Breeder, the Thick-Billed Murre, at a Breeding Colony in the Canadian Arctic[END_REF][START_REF] Albert | Using Blood and Feathers to Investigate Large-Scale Hg Contamination in Arctic Seabirds: A Review[END_REF] and differences in regional Hg contamination.

Environmental variables explaining Hg concentration in kittiwakes

Among the environmental variables selected to explain Hg concentrations in kittiwakes, Chl a concentration was the strongest predictor. Extreme Chl a concentrations (i.e. lowest and highest events) were associated with higher Hg exposure in kittiwakes. Chl a concentrations remain a reliable predictor of phytoplankton biomass and thus primary productivity in open-water [START_REF] Huot | Does Chlorophyll a Provide the Best Index of Phytoplankton Biomass for Primary Productivity Studies?[END_REF] . New primary production is directly available to mesozooplankton which in turn supports production of upper trophic levels 36 . This suggests that extremely high productivity and extremely low productivity would result in high Hg concentrations in kittiwakes. Interestingly, Chl a and water temperature are strong predictors of shifts in fish communities [START_REF] Eisner | Pelagic Fish and Zooplankton Species Assemblages in Relation to Water Mass Characteristics in the Northern Bering and Southeast Chukchi Seas[END_REF][START_REF] Mérillet | Environment Outweighs the Effects of Fishing in Regulating Demersal Community Structure in an Exploited Marine Ecosystem[END_REF][START_REF] Collie | Long-Term Shifts in the Species Composition of a Coastal Fish Community[END_REF] . For example, in the Bering and Chukchi seas, the distributions and community composition of pelagic species were dependent of temperature and chlorophyll concentrations [START_REF] Eisner | Pelagic Fish and Zooplankton Species Assemblages in Relation to Water Mass Characteristics in the Northern Bering and Southeast Chukchi Seas[END_REF] . Polar cod Boreogadus saida were most abundant in warm Alaska Coastal Water of the central Chukchi Sea, characterized by low chlorophyll, low nutrients, and small zooplankton taxa. In contrast, adult Pacific herring Clupea pallasii were more abundant in the less stratified Bering Strait waters and in the colder Bering Shelf Water of the northern Bering and southern Chukchi seas, characterized by high chlorophyll, high nutrients, and larger zooplankton taxa. Among the five selected environment-related variables, Chl a was the best predictor of trophic level and carbon source in kittiwakes' blood, this confirms the fact Chl a could be a reliable indicator of fish community shifts.

Kittiwakes are generalist species that can shift trophic level according to environmental conditions [START_REF] Pratte | Variable Sea-Ice Conditions Influence Trophic Dynamics in an Arctic Community of Marine Top Predators[END_REF] and a shift in trophic level has been related to a modification of Hg exposure in other Arctic seabirds Uria lomvia 6 . As previously mentioned, since 2005, Kongsfjorden kittiwakes have been shifting their polar cod dominated diet by a more generalist diet composed of species from sub-Arctic and Atlantic origin, such as capelin (Mallotus villosus), Atlantic herring (Clupea harengus), young Atlantic cod (Gadus morhua) and haddock (Melanogrammus aeglefinus) 31 . Sub-Arctic and Atlantic pelagic fish have higher temperature tolerances and generally benefit from warmer water temperatures [START_REF] Kunz | New Encounters in Arctic Waters: A Comparison of Metabolism and Performance of Polar Cod (Boreogadus Saida) and Atlantic Cod (Gadus Morhua) under Ocean Acidification and Warming[END_REF] . Interestingly, Hg concentrations in polar cod muscles sampled in the East of the Barents Sea were 2.8 fold higher than in capelin, whereas in Atlantic cod, Hg concentrations were 3.2 fold higher than in capelin [START_REF] Gopakumar | Mercury in Barents Sea Fish in the Arctic Polar Night: Species and Spatial Comparison[END_REF] . Similar results were reported in the Canadian Arctic, with total Hg in polar cods that were higher than in capelin [START_REF] Pedro | Mercury and Persistent Organic Pollutants in Native and Invading Forage Species of the Canadian Arctic: Consequences for Food Web Dynamics[END_REF] .

From 2013 to 2019, the occurrence of Atlantic herring in Kongsfjorden kittiwake's diet has often been above 20% (Eulaers I., Gabrielsen G.W. et al., personal communication), and Hg concentrations are similar between polar cod and herring [START_REF] Jaeger | Biomagnification of Mercury in Selected Species from an Arctic Marine Food Web in Svalbard[END_REF] . Consequently, a progressive replacement of polar cod to capelin, then capelin to other more Hg contaminated species such as Atlantic cod juveniles or Atlantic herrings in Kongsfjorden kittiwakes diet could partly explain the U-shaped relationship in kittiwakes' blood Hg concentrations from 2000 to 2019.

Sea-ice, which was the strongest predictor of non-detrended Hg concentrations in kittiwakes, plays a large role in the accumulation of GEM, and acts as a barrier for GEM air-sea exchange, thus it has been proposed that melting sea-ice in a warming climate would lead to increased GEM evasion, resulting in decreased Hg levels in the surface Arctic Ocean [START_REF] Dimento | The Impact of Sea Ice on the Air-Sea Exchange of Mercury in the Arctic Ocean[END_REF] . In a recent study, Schartup and colleagues [START_REF] Schartup | Influence of the Arctic Sea-Ice Regime Shift on Sea-Ice Methylated Mercury Trends[END_REF] observed that first-year ice contains more MeHg than multi-year ice and calculated that the ongoing shift from an old multi-year ice to a first-year ice dominated Arctic Ocean may cancel out the previously expected decline in MeHg contributions from seaice [START_REF] Schartup | Influence of the Arctic Sea-Ice Regime Shift on Sea-Ice Methylated Mercury Trends[END_REF] . Therefore, although a reduction in sea-ice extent would accelerate GEM evasion and sunlight-induced breakdown of MeHg in surface waters [START_REF] Point | Methylmercury Photodegradation Influenced by Sea-Ice Cover in Arctic Marine Ecosystems[END_REF] , these mechanisms could be counterbalanced by the release of MeHg from young sea-ice.

In a study on two other seabird species, the thick-billed murre and the northern fulmar in the Canadian High Arctic, a 40-year time series showed that NAO and surface air temperature were the strongest predictors of Hg concentrations in eggs 21 . Snowfall and sea-ice concentrations were also strong predictors for murres and fulmars, respectively, however these relationships were only observed after time lags between four (NAO, summer temperature) and seven years (sea-ice concentrations). It is interesting to see that, in Svalbard, environmental conditions in the preceding spring can immediately influence blood Hg concentrations in kittiwakes. This could result from the labile characteristic of blood and the numerous biogeochemical exchanges that occur on the western coast of Spitsbergen. The area around Spitsbergen and Fram Strait is one of the most climatically sensitive areas in the world 3,[START_REF] Rogers | The Role of Fram Strait Winter Cyclones on Sea Ice Flux and on Spitsbergen Air Temperatures[END_REF] . Therefore, in comparison to the high Canadian Arctic, Hg exchanges between the environment and biota could be more dynamic and, thus, accelerated.

We have actually no evidence that those increased concentrations of Hg in kittiwake blood will result in deleterious effects on colonies' size. In more recent years some colonies seem to have recovered from the decline observed between 1970 and 2015 [START_REF] Descamps | Species and Spatial Variation in the Effects of Sea Ice on Arctic Seabird Populations[END_REF] . Atlantification or borealization processes in some feeding habitats may be offering a rather favourable environment to kittiwakes [START_REF] Descamps | Species and Spatial Variation in the Effects of Sea Ice on Arctic Seabird Populations[END_REF] and the toxic effects of Hg could be counterbalanced by a higher nutritive value of Atlantic prey.

To conclude, from 2000 to 2013, blood Hg in breeding kittiwakes from a Svalbard colony decreased by 3% per year, and from 2014 to 2019, blood Hg increased by 11% per year. Sex was the strongest biological predictor of Hg concentrations in kittiwakes, whereas the strongest environmental predictor of Hg in kittiwakes was Chl a, with extreme values associated to high Hg concentrations. Chl a, but also sea-ice extent are reliable predictors of fish communities.

The blood Hg U-shaped trend is likely due to a shift from polar cod to Atlantic fish species in kittiwake diet, likely added in the later years to the input of MeHg from Arctic shrinking seaice. Moreover, our results raise concern for Hg contamination in the wintering ground of kittiwakes, because the concentrations measured in pre-laying birds have been associated to decreased breeding probability and disrupted parental care in this species.

Table 1. Biological factors influencing blood Hg concentrations in black-legged kittiwakes from Kongsfjorden Svalbard.

A) Blood Hg concentration (n=1066) and stable isotope values (n=80) across the breeding cycle. B) Blood Hg variations according to biometrics (n=1066). Models are LMMs computed on detrended variables, with year and bird identity as random factors. Values are estimates and 95% confidence intervals. The subset of recaptured birds consists of 43 samples representing 20 individuals recaptured within the same year between two and three times. In birds, MeHg storage from blood into tissues (e.g. muscles) represents 50-95% of the Hg ingested for up to a period of approximately 8 ± 2 weeks 1,2 . When birds arrive to the breeding ground, their blood Hg concentrations therefore reflect the exposure of the wintering grounds 3 . The eight weeks prior to the chick-rearing period include the pre-laying period (the courtship and nest building stage), where females forage exclusively inside the fjord, whereas males forage inside and outside the fjord at the Greenland-Svalbard ridge and the Fram Strait 4 . In bird red blood cells, the turnover rate of δ 13 C and δ 15 N values (reflecting feeding habitat and trophic level, respectively) is around 30-50 days 5,6 , thus sampling conducted during the chickrearing stage represents the feeding habits between the pre-laying and incubating periods.

A) Hg and stable isotopes variation according to breeding phase

Molecular sexing, mercury and stable isotope analyses

According to the planned subsequent analyses, blood samples were centrifuged or not, therefore samples comprised either whole blood or red blood cells, which were stored at -20°C until analyses. Molecular sexing was performed at the Centre d'Etudes Biologiques de Chizé (CEBC) on red blood cells as detailed elsewhere 7 . Prior to Hg and isotopic analyses, red blood cells and whole blood samples (hereafter called blood) were freeze-dried and powdered. We measured total Hg (THg) in 1049 kittiwakes, at Littoral Environnement et Sociétés lab (LIENSs) as previously described 8 . MeHg represents more than 90% of THg in seabirds 9 . We used an Advanced Mercury Analyzer spectrophotometer (Altec AMA 254) with a detection limit of 0.1 ng Hg. At least two aliquots ranging from 5 to 10 mg were analyzed for each samples and quality assessment was measured by repeated analyses (every 10 samples) of certified reference material (CRM) TORT-2 (lobster hepatopancreas, NRCC; certified value 0.27±0.06 µg/g). CRM measurements were in good agreement with the certified values and recoveries were 103.7±0.01%. Hg concentrations are expressed in µg/g dry weight (dw). Hg concentrations obtained from whole blood and red blood cells were close to identical (R 2 =0.98).

Stable isotopes were analysed in 405 kittiwakes. The isotopic niche of kittiwakes was used as a proxy of their trophic niche, as blood cells is a metabolically active tissue that integrates periods of weeks before sampling 5,10 . The isotopic method was validated in seabirds, with δ 13 C values indicating their foraging habitats 11,12 and their δ 15 N values increasing with trophic levels 13 . Sub-samples of homogenized blood powder were then weighed (~0. The consistently low C:N values of blood cells verified that low lipid content did not necessitate lipid extraction 14 . The isotopic ratio of red blood cells slightly differed from that of whole blood, we therefore applied a correction factor to stable isotope values obtained from whole blood samples (red blood cell δ 13 C =whole blood δ 13 C × 0.985 and red blood cell δ 15 N = whole blood δ 15 N × 1.004).

Ocean acidification indicators

Aragonite saturation state is commonly used to track ocean acidification because it is a measure of carbonate ion concentration ([CO3 2-]). At a lower pH, metals that form strong complexes with hydroxide (OH -) and CO3 2-will be more abundant in their free forms (e.g. Hg0). These changes in speciation also increase the thermodynamic and kinetic activity of the metals (Millero, 2001a; b). As the pH of seawater decreases, metals, such as Hg that forms strong complexes with OH -and CO3 2-, will undergo significant changes in speciation [15][16][17][18] . Adsorption of metals into organic material will also change at a lower pH. Ocean acidification measures (pH and aragonite saturation) were measured in July in Kongsfjorden and were not available from 2000-2012 19 . Since pH and aragonite saturation were measured in July and the Hg concentrations in blood referred to bird exposure since May, we associated ocean acidification measures from year n-1 to blood Hg concentrations from year n, to avoid testing predictors obtained 'a posteriori' on our response variable.

Statistical analyses

First, we investigated which biological factors (i.e. sex, body-size/-condition, or feeding habits)

would influence individual Hg concentrations. To do so, we analysed the variations of Hg across the breeding cycle (pre-laying, incubation, and chick-rearing periods) and its interaction with sex, by using linear mixed-effect models (LMMs) with breeding period and sex as independent categorical 3-and 2-level predictors and bird identity as a random factor to explain detrended blood Hg concentrations. We obtained estimated marginal means (EMMs) for specified factor combinations of the LMM (breeding period × sex interactions) by using the 'emmeans' package 20 . To confirm the obtained trend, we analysed repeated sampling of the same individuals across the breeding cycle in a subset of the largest dataset, by including only individuals that had been captured at least twice in the same year. In addition, we analysed changes in feeding habits across the breeding cycle. To do so, we used LMMs with δ 13 C and δ 15 N (proxies of feeding habitat and trophic level, respectively) as response variables, breeding period as a categorical 3-level predictor, and bird identity as a random factor. Sufficient data for each period of the breeding cycle were only available in 2011.

Then, we tested whether body measurements, could predict Hg concentrations in kittiwakes, using LMMs with bird identity as a random factor. Mercury concentration was set as the response variable and body mass, size corrected mass SCM (i.e. residuals from the mass ~ skull size regression), scaled mass index SMI 21 and skull size were set as predictors in different candidate models because most predictors were correlated. We ranked the candidate models using Akaike Information criterion (AIC) 22 , the best model (lowest AIC) was considered as the most reliable biological predictor of detrended Hg concentrations.

Second, to investigated the temporal trends of Hg and other biological variables in chick-rearing kittiwakes, over a 20-year period, we used GAMMs on raw data to account for non-linear relationships between response variables and predictors, and for recaptured individuals. Year, δ 13 C or δ 15 N and sex were used as continuous predictors, and bird identity as a random factor to analyse temporal trends of Hg. Since Hg concentrations presented non-linear variations according to year, we entered year as a smoothed term (k=4) using penalized regression splines and a Gaussian distribution 23 . We used an identical approach (GAMM) to study the temporal trends of δ 13 C, δ 15 N, body mass and skull size, with year and sex as continuous predictors, and bird identity as a random factor. We adapted the smoothed term to avoid overfitting. We ranked the candidate models using Akaike Information criterion (AIC) 22 , the best model (lowest AIC) was considered as the most reliable biological predictor of raw Hg concentrations. Break-points were identified by using a Bayesian approach ('broken stick' modelling technique 24,25 , based on beta regression analysis with 'betareg' package in R 26 . This method infers the 'year' where blood Hg concentration starts to translate statistically as a change-point in the distribution of blood Hg concentration versus year. To estimate the year of the change-point Ti, in blood Hg concentrations, we used a profile likelihood approach where the likelihood was computed for each year. The value of Ti that maximized the likelihood was thus evaluated, and an approximate confidence interval for Ti was computed with a Likelihood Ratio Test with 1 df.

To avoid pseudoreplication due to the large, and sometimes unbalanced, number of kittiwakes caught per year, for each year, we calculated the average Hg concentrations, δ 13 C and δ 15 N values for each sex, and performed statistical analyses on this subset of n=38. We tested both linear (e.g. sea ice extent) and quadratic (e.g. sea ice extent + sea ice extent 2 ) relationships using

  (aragonite: 0.90 (-0.96; 2.76), pH: -2.34 (-11.8; 7.10), values are estimates and 95% CIs obtained from LMMs).

Figure 1 .

 1 Figure 1. Blood Hg concentrations and stable isotope values across the breeding cycle in black-legged kittiwakes (Kongsfjorden Svalbard). A) Boxplot of non-adjusted Hg concentrations (n=1049, 2000-2016) in males (blue) and females (pink) for each breeding period (PL-Pre-laying, INC-Incubation, CR-Chick-rearing). Median is shown with lower and upper quartiles (box), whiskers (1.5 × interquartile range from box) and outliers (values beyond whiskers), B) Mean nitrogen (δ 15 N) and carbon (δ 13 C) stable isotope values with error bars (SE) in blood samples from black-legged kittiwakes (n=77, 2011) during each breeding period.

Figure 2 .

 2 Figure 2. Temporal trends (2000-2019) of blood Hg concentrations in adult chick-rearing kittiwakes (n=392) from Kongsfjorden Svalbard. Hg partial residuals are derived from the best model selected by AIC (Table 2), Hg partial residuals are adjusted for δ 15 N, sex as fixed factors, and bird identity as a random factor. The red line represents the selected change-point in Hg concentrations which corresponds to 2014. Values are estimates and 95% CI obtained from LMMs for the period of 2000-2013 on the left and 2014-2019 on the right.

Figure 3 .

 3 Figure 3. Temporal trends (2000-2019) of stable isotope values in red blood cells and biometrics in adult chick-rearing kittiwakes (n=392) from Kongsfjorden Svalbard. A-C)Feeding habits are inferred from nitrogen and carbon stable isotope values in blood (δ15 N and δ 13 C, respectively), biometrics include B) body mass and D) skull length. The red line represents the year selected as a change-point for temporal variations of Hg in blood of kittiwakes. Dots are partial residuals derived from a GAMM including year as a smoother, sex as a fixed factor, and bird identity as a random factor. Diagnostic plots are given in supplementary information (FigureS1).

Figure 4 .

 4 Figure 4. Influence of chlorophyll a on blood Hg concentrations and feeding habits in chick-rearing kittiwakes between 2000-2019, (n=38, Kongsfjorden, Svalbard). Feeding habits are inferred from nitrogen and carbon stable isotope values in blood (δ 15 N and δ 13 C, respectively). Dots are year average values for each sex (females in red and males in blue), the lines are the quadratic and linear fits to the predicted values and the shaded areas represent their 95% confidence interval derived from LMs on detrended variables. Females are in red and males in blue.

  4 mg) with a microbalance and packed in tin containers. Relative abundances of C and N isotopes were determined with a continuous flow mass spectrometer (Thermo Scientific Delta V Advantage) coupled to an elemental analyzer (Thermo Scientific Flash EA 1112). Isotopic results are presented in the δ* notation relative to Vienna PeeDee Belemnite and atmospheric N2 for δ 13 C and δ 15 N, respectively. Replicate measurements of internal laboratory standards (acetanilide) indicate measurement errors < 0.10 ‰ for both δ 13 C and δ 15 N values. The C:N mass ratios of the samples were calculated as the ratio between the mass percentages in carbon and nitrogen.

Figure S1 .

 S1 Figure S1. Temporal variations of climate-related variables measured in spring in Svalbard or the Fram Strait. The blue line represents a quadratic relationship between the environmental variables and year, and the grey shaded area the standard error.

Figure S2 .

 S2 Figure S2. Blood Hg variations across the breeding cycle in kittiwakes captured several time. Each dot connected by a line represents Hg variation across the cycle, males are in green and females in pink (LMM: intercept 1.97 (1.59; 2.35), incubating birds -0.43 (-0.77; -0.10), p=0.021, chick-rearing -0.81 (-1.12; -0.50), p<0.001).

Figure S3 .

 S3 Figure S3. Diagnostic plots for the five GAM models investigating the biological variables which better explain temporal trends of Hg in kittiwakes. Diagnostic plots determine the goodness of the fit of the model to the data. For each model is produced a quantile-quantile plot of residuals (Q-Q plot), a histogram of residuals, a 'residuals vs linear predictor' plot, and an 'observed vs fitted values' plot.

Figure S4 .

 S4 Figure S4. Diagnostic plots for GAM models investigating the temporal trends of stable isotopes and biometrics in kittiwakes. Diagnostic plots determine the goodness of the fit of the model to the data. For each model is produced a quantile-quantile plot of residuals (Q-Q plot), a histogram of residuals, a 'residuals vs linear predictor' plot, and an 'observed vs fitted values' plot.

Table 2 . Predictors of blood Hg variation in chick-rearing kittiwakes (n=392) from Kongsfjorden, Svalbard (2000-2019). Model selection and model-averaged estimates.

 2 The first part of the table shows the different candidate models and their AIC ranking, the most competitive model appears on the first line (ΔAIC = 0). The second part of the table shows the estimates and 95% confidence intervals derived from GAMMs. Bold values represent significant relationships. GAMMs including year as a smoother, sex, body mass, skull size, δ 15 N, δ13 C and sex as fixed factors and bird identity as a random factor. Sex was included in all models except for those with body measurements which are already sex-dependent.

	AIC selection
	Intercept	s(Year, k=5) Sex Mass	Skull	δ 15 N δ 13 C df logLik AIC ΔAIC weight
	-0.142	+	+		0.09078	7 -142.955 299.9	0	0.817
	3.447	+	+		0.1193 7 -144.661 303.3 3.41 0.148
	0.9943	+	+			6 -147.175 306.3 6.44 0.033
	-4.252	+		0.04553 0.09819	7 -148.938 311.9 11.97 0.002
	-0.6652	+		0.04935	0.1312 7 -150.556 315.1 15.2	0
	-3.364	+		0.0494		6 -153.966 319.9 20.02	0
	-1.351	+	0.002474		0.1232	7 -163.333 340.7 40.76	0
	2.537	+	0.002289		0.11 7 -170.93 355.9 55.95	0
	0.2107	+	0.002465			6 -172.255 356.5 56.6	0
	1.154	+				5 -174.943 359.9 59.98	0
	1.164					3 -199.323 404.6 104.74	0
	Conditional model averaged estimates
	Variable				β (95% CI)
	Intercept				0.42 [-2.24; 3.07]
	δ 15 N				0.09 [0.04; 0.14]
					0.29 [0.22; 0.37]
					0.12 [0.04; 0.19]
					0.05 [0.03; 0.06]
				0.002 [0.001; 0.003]

Sex

δ

13

C

Skull

Mass

The '+' sign in the upper part of the table indicates that a smoothed term "s(Year, k=5)" or categorical term "Sex" was included in the model.

Table S1 . Overview of the samples available per year and per breeding stages.

 S1 Values are means and standard deviation for blood Hg concentrations (μg/g dw), nitrogen δ 15 N and carbon δ13 C stable isotope values (‰), body mass (g) and skull length (mm). For biometrics, in some years the number of kittiwakes weighed differed from that of measured, in that case the first number (N1) in (N1/N2) from the biometrics column represents the number of kittiwakes weighed and N2 the number of kittiwakes measured. Non-available data = na. All samples from year 2010 were missing.

			Mercury (μg/g dw)	δ	15 N (‰) δ	13 C (‰)	Biometrics		Body mass (g)	Skull length (mm)
				Both sex			Both sex	Females Males	Females	Males	Females	Males
	Year	Breeding period N Mean SD	N Mean SD Mean SD	N	N Mean SD Mean SD Mean SD Mean SD
	2000	Chick-rearing	17 1.51 0.50	17 13.0 0.4 -20.0 0.2	10	7	367 27.9 403 11.5 87.4 1.9 93.0 2.7
	2001	Chick-rearing	30 1.59 0.32	19 13.1 0.5 -20.0 0.2	13/15 13/15 347 21.3 384 32.9 88.4 2.6 94.4 1.7
	2002	Chick-rearing	28 1.24 0.37	20 12.4 0.6 -20.6 0.2	15/10 13/9 368 25.6 408 29.8 88.7 1.3 93.4 1.5
	2003	Chick-rearing	20 0.78 0.25	20 11.9 0.6 -21.5 0.4	10	10	357 18.9 379 44.2 88.9 1.1 92.7 1.9
	2004	Chick-rearing	30 1.44 0.37	20 12.8 0.3 -20.5 0.1	15	15/14 392 24.2 425 22.2 89.5 1.5 94.0 2.1
	2005	Chick-rearing	31 0.99 0.24	7 14.0 0.3 -20.2 0.1	11/7 20/16 353 15.1 407 22.4 90.6 1.6 94.2 1.5
	2006	Chick-rearing	29 1.38 0.42	20 12.9 0.6 -20.7 0.3	13	16/15 347 29.9 383 17.4 88.5 1.1 93.6 1.5
	2007	Chick-rearing	30 0.86 0.21	20 12.7 0.3 -20.0 0.2	15/13 15/13 381 28.8 432 25.9 89.4 1.9 93.2 1.8
	2008	Pre-laying	69 2.03 0.44	0 na na na na	31/28 38/34 396 31.7 423 32.1 88.9 1.2 93.6 2.3
	2008	Incubation	38 1.54 0.46	20 11.8 0.5 -20.3 0.1	15	23/22 395 16.4 440 27.7 89.7 1.5 93.3 1.9
	2008	Chick-rearing	34 0.90 0.31	16 12.1 0.4 -20.3 0.2	20/18 14/13 354 16.8 411 22.4 89.2 1.5 93.0 1.9
	2009	Pre-laying	46 2.18 0.51	0 na na na na	20/15 26/21 388 36.6 427 31.6 89.2 1.9 94.0 1.9
	2009	Incubation	2 2.00 0.60	0 na na na na	na	na	na	na	na	na	na na na na
	2009	Chick-rearing	17 0.96 0.45	10 12.5 0.6 -20.9 0.5	9	8	346 17.5 385 24.8 89.3 2.0 94.2 1.9
	2010	NA	0	na	na	0 na na na na	na	na	na	na	na	na	na na na na
	2011	Pre-laying	93 1.85 0.46	30 14.1 0.3 -20.2 0.3	39	53/54 394 36.3 413 25.6 88.8 1.8 93.9 1.7
	2011	Incubation	52 1.58 0.62	28 13.0 0.4 -20.4 0.1	24	29	382 25.6 423 28.1 89.1 1.8 93.0 2.1
	2011	Chick-rearing	34 0.87 0.24	22 12.9 0.6 -20.6 0.2	18	16	338 18.2 376 24.8 88.5 1.7 92.8 1.9
	2012	Pre-laying	3 1.80 0.30	0 na na na na	na	na	na	na	na	na	na na na na
	2012	Incubation	102 1.52 0.51	0 na na na na	51	51	384 22.3 429 28.0 88.0 1.9 92.9 1.8
	2012	Chick-rearing	51 1.09 0.37	43 13.4 0.5 -20.0 0.2	26	25	348 25.8 394 22.8 88.1 1.4 93.1 1.5
	2013	Incubation	80 1.39 0.34	2 11.5 0.9 -20.1 0.2	23	57	385 20.3 429 25.2 88.7 1.3 93.7 1.7
	2013	Chick-rearing	53 1.22 0.33	18 11.3 0.4 -19.9 0.1	11	42/41 374 27.3 414 27.2 87.3 1.6 93.7 1.9
	2014	Incubation	4 1.49 0.43	0 na na na na	2	2	388 10.6 453 38.9 87.9 0.1 94.5 0.1
	2014	Chick-rearing	34 0.77 0.31	20 11.9 0.4 -20.0 0.1	15	19	402 19.5 451 27.1 89.4 2.1 93.7 1.6
	2015	Incubation	68 1.69 0.58	0 na na na na	34	34	378 25.2 409 17.2 89.0 1.7 93.8 1.8
	2015	Chick-rearing	31 0.88 0.33	31 11.7 0.3 -21.3 0.2	16	15	371 19.7 402 31.0 90.0 2.4 94.3 1.7
	2016	Pre-laying	0	na	na	0 na na na na	na	50	na	na 417 22.8 na na 93.9 1.6
	2016	Incubation	0	na	na	0 na na na na	24	18	380 24.2 412 18.5 88.6 2.2 93.2 2.1
	2016	Chick-rearing	22 1.04 0.29	22 12.4 0.6 -20.6 0.2	11	11	360 25.4 409 21.3 90.7 1.2 95.1 1.8
	2017	Chick-rearing	33 1.45 0.32	33 12.7 0.6 -21.1 0.2	22	15	378 27.6 433 31.5 90.6 1.2 95.3 1.7
	2018	Chick-rearing	46 1.48 0.35	46 12.9 0.1 -20.6 0.1	23/22	24	351 20.3 386 18.6 89.1 1.4 93.3 1.5
	2019	Chick-rearing	24 1.16 0.27	24 13.2 0.5 -21.0 0.1	10	9	353 13.2 394 27.4 89.9 1.0 94.1 2.0

Table S5 . Effects of environmental variables on Hg concentrations in Svalbard chick- rearing kittiwakes (n=38).

 S5 The table shows the ten most competitive linear models (ranked using AIC) among a list of 34 for both raw (upper part) and detrended variables (lower part).

	Model Intercept Sex	poly (Sea Ice, 2)	poly (Chl a, 2)	poly (SSTA, 2)	poly (Air Temp, 2)	poly (Precip, 2)	poly (GEM, 2) NAO	Sea Ice	Chl a SST A Precip GEM df logLik AIC ΔAIC weight
	Raw variables
	13	0.9776	+	+						5 4.99	0	0	0.995
	17	0.9776	+					+		5 -1.89 13.8 13.8 0.001
	30	0.6937	+							0.3542	4 -3.48 15 15.0 0.001
	5	1.139		+						4 -3.65 15.3 15.3	0
	31	0.9621	+							0.3677	4 -3.76 15.5 15.5	0
	15	0.9776	+			+				5 -2.96 15.9 15.9	0
	33	0.8306	+							0.00166	4 -4.02 16 16.0	0
	34	1.906	+							-0.6416 4 -4.16 16.3 16.3	0
	14	0.9776	+		+					5 -3.27 16.5 16.5	0
	1	0.9776	+							3 -5.42 16.8 16.8	0
	Detrended variables
	14	-0.1884	+		+					5 7.14 -4.3	0	0.231
	31	-0.1828	+							0.2643	4 5.25 -2.5 1.8 0.094
	1	-0.1884	+							3 4.16 -2.3 2.0 0.086
	34	-0.1987	+							-0.6063 4 5.03 -2.1 2.2 0.076
	33	-0.1802	+							0.00101	4 4.96 -1.9 2.4 0.071
	18	-0.1884	+						+	5 5.74 -1.5 2.8 0.057
	27	-0.186	+						0.07604	4 4.73 -1.5 2.8 0.056
	16	-0.1884	+				+			5 5.52	-1	3.3 0.045
	13	-0.1884	+	+						5 5.43 -0.9 3.4 0.042
	29	-0.1875	+							-6E-04	4 4.36 -0.7 3.6 0.039

ACKNOWLEDGEMENTS

This project was supported by Institut Polaire Français (project n° 330, O. Chastel), the ANR ILETOP (ANR-16-CE34-0005) and the program SENSEI funded by the BNP Paribas Foundation (C. Barbraud). Additional funding was provided by the Fram center for fieldwork (Multiple Stressor seabird project). The authors thank the CPER and the FEDER for funding the AMA and IR-MS devices; K. Sagerup, S. Nilsen, A. Ask, T. Nordstad, and C. Clément-Chastel for their help in fieldwork, B. Picard for his help with the database, C. Churlaud and M. Brault-Favrou for their assistance in Hg analyses, G. Guillou for running isotope analyses and H. Routti for her help in statistical analyses. The authors acknowledge K. Aspmo

Pfaffhuber, NILU, EMEP and Norwegian Environment Agency for providing data and funding these measurements. The IUF (Institut Universitaire de France) is acknowledged for its support to P. Bustamante as a Senior Member. The authors thank S. Hope and three anonymous reviewers for improving the English and the manuscript. We declare no conflict of interest.

The data that support the findings of this study are available from the corresponding author upon reasonable request.

Table 3. Relationship between detrended blood Hg concentrations in chick-rearing kittiwakes and environmental factors (n=38).

The table shows the 5 highest AIC ranked models, the most competitive model appears on the first line (ΔAIC = 0). We give the model averaged estimates and 95% confidence intervals derived from LMs with Hg concentrations as a response variable. Bold values represent significant relationships. The complete list of candidate models and results from non-detrended variables are given in the supplementary information.

Explanatory variables df logLik AICc ΔAIC weight Linear Quadratic

Chl a + Chl a 2 + Sex 5 7. 14 Chl a : Chlorophyll a , GEM: Gaseous Elemental Mercury, Air Temp: Air temperature, Sea Ice: Sea ice extent, SSTa: Sea surface temperature anomaly.

AIC ranking

Supporting information MATERIAL AND METHODS

Fieldwork and blood sampling

From 2000 to 2016 (excluding 2010), 1049 kittiwakes have been captured during different phases of the breeding cycle at the Krykkjefjellet "main colony" in Kongsfjorden, Svalbard (78°54'N, 12°13'E), 7 km southeast of Ny-Ålesund, Norway. We studied kittiwakes in one plot of around 117 pairs breeding on cliff ledges at heights of 5-10 m. Kittiwakes were caught on their nests with a noose at the end of a 5 m fishing rod. A blood sample (ca. 1 mL) was collected from the alar vein with a 1 mL heparinised syringe and a 25-gauge needle. Kittiwakes were individually marked with metal rings and PVC plastic bands engraved with a three-digit code and fixed to the bird's tarsus for identification from a distance. Birds were weighed to the nearest 5 g using a Pesola spring balance, and their skull length was measured to the nearest 0.5 mm with a sliding calliper, biometrics were available for 997 kittiwakes. In several years, birds were captured across the breeding period, which can be separated into three phases: the pre-laying period from end-April to early June which corresponds to the arrival to the colony, nest building and mating; the incubation period lasts from early/mid-June to mid-July and starts when the first egg is laid, then comes the chick-rearing period (~mid-July to mid-August) which starts when the first egg has hatched and lasts until all the chicks have fledged. In several occasions, within a year, individuals were captured at least during two stages (n=17), and a few times during the three stages (n=3), enabling us to study blood Hg dynamic across the breeding cycle.

Temporal trends of blood Hg have been analysed in chick-rearing kittiwakes only, in order to reflect local Hg exposure and to make between year comparisons possible.

Mercury and stable isotopes in blood

linear models (LMs), since Hg trends in biota often present non-linear variations [40][START_REF] Thiebot | Gentlemen First? 'Broken Stick' Modelling Reveals Sex-Related Homing Decision Date in Migrating Seabirds[END_REF][START_REF] Gaden | Mercury Trends in Ringed Seals (Phoca Hispida) from the Western Canadian Arctic since 1973: Associations with Length of Ice-Free Season[END_REF] and have been increasing non-linearly in Svalbard vertebrates 17,18 . We computed a list of 34 biologically relevant model candidates by adding the variables selected from the best model from the analysis described on step two (Hg variations in relation to morphometrics and sex), to which we added one of the eight climate-related variables without missing years, either as linear or quadratic. We show the ten most competitive models in Table S3 for both detrended and raw data. We then selected the best fits (linear or quadratic) for all variables included in the top 10 models to explain detrended Hg concentrations, and used conditional model averaging to make inference from all of the models [START_REF] Rigét | Temporal Trends of Hg in Arctic Biota, an Update[END_REF][START_REF] Zabala | Early Breeding Failure in Birds Due to Environmental Toxins: A Potentially Powerful but Hidden Effect of Contamination[END_REF] . Variables with missing years (pH and aragonite saturation at the bottom of the water column) were analysed separately using LMs. Finally, we tested which of the strongest climate-related variables would better explain feeding habits (δ 13 C or δ 15 N), to do so we used LMs with detrended variables.

All models were ranked using AIC [START_REF] Burnham | Model Selection and Multimodel Inference: A Practical Information-Theoretic Approach[END_REF] . We obtained the number of parameters (K), the difference in AIC values between the "best" model and the model at hand (ΔAIC), and a normalized weight of evidence in favor of the specific model, relative to the whole set of candidate models, derived by e (-0.5(ΔAICc)) (AICc weights). We obtained estimates and 95% confidence intervals (95% CI) for each parameter, and an effect was considered significant when the 95% CI did not include 0. Model fit was assessed by using residual diagnostic plots. All analyses and figures were performed using R software v.3.6.3 [START_REF] Lukacs | Model Selection Bias and Freedman's Paradox[END_REF] with the following packages: car, mgcv, nlme, effects, ggplot2 and corrplot [START_REF] Akaike | Factor Analysis and AIC[END_REF][START_REF] Core | R: A Language and Environment for Statistical Computing[END_REF][START_REF] Fox | An {R} Companion to Applied Regression[END_REF][START_REF] Fox | Visualizing Fit and Lack of Fit in Complex Regression Models with Predictor Effect Plots and Partial Residuals[END_REF][START_REF] Pinheiro | Nlme: Linear and Nonlinear Mixed Effects Models[END_REF][START_REF] Wei | Corrplot": Visualization of a Correlation Matrix[END_REF] .

Data sharing policy

The data that support the findings of this study are available from the corresponding author upon reasonable request.