Introduction

The freshwater cycle in high-latitude seas is highly seasonal. Most of the river runoff, precipitation, and sea ice melt is added in spring and summer, while in winter, river runoff typically is low as the land remains frozen, and the formation of sea ice effectively withdraws freshwater from the surface ocean. This seasonal freshwater cycle exerts a first-order control on the oceanography of the Arctic’s seas (Carmack and Macdonald 2002; Carmack et al. 2016). In summer, the addition of runoff, precipitation, and sea ice melt lead to strong stratification. In winter, the formation of sea ice adds cold and salty brine to the surface water, which reduces stratification and promotes vertical mixing (Granskog et al. 2011). However, freshwater remaining at the end of summer and/or addition of freshwater runoff during winter may counter-act the addition of brine from sea ice growth and thus maintain shallow stratification. The balance between freshwater supply and the production of brine by sea ice therefore controls the depth and intensity of ocean convection (e.g., Macdonald 2000). The central role of freshwater seasonality at local, regional, and system-wide scales raises concern about the oceanographic consequences of shifting river inflow from summer to winter, as predicted for Arctic seas under a warmer climate, and through hydroelectric regulation, which, in some cases, stores water in spring and releases it in winter (Déry et al. 2011; Prinsenberg 1991).

Hudson Bay in northern Canada is a large semi-enclosed shelf sea (Fig. 1) similar to the Arctic Ocean in receiving large freshwater inputs from river runoff and the seasonal melt of sea ice. The Bay is undergoing very rapid change in its freshwater balance. First, Hudson Bay is experiencing climate change as an earlier break up and a later freeze up in its seasonal sea ice cycle (Hochheim and Barber 2014). The trend towards a longer open-water season is expected to continue into the future (Gagnon and Gough 2005a; Huard et al. 2014). Second, total annual river discharge to Hudson Bay has increased since 1990 (Déry et al. 2011, 2016) and an additional, more rapid, increase in discharge is projected for the future particularly in fall and winter (Collins et al. 2013; Gagnon and Gough 2005b). Third, there has been a shift in streamflow towards higher flows in winter (see inset in Fig. 1) particularly in eastern Hudson Bay and James Bay (Déry et al. 2011, 2016). The increase in winter river discharge is partly due to the development after the late-1970s of the James Bay hydroelectric complex and the La Grande system, which includes several large reservoirs where water can be retained in spring and summer and released in winter for hydropower production to meet peak demands in central Canada (Déry et al. 2016).

Fig. 1
figure 1

Map of Hudson and James Bay showing the general cyclonic circulation pattern of surface waters and location of the Belcher Islands study area (red star is Sanikiluaq). The upper inset shows changes in the seasonal hydrograph for the whole Hudson and James Bay system (23 gauged rivers) between the periods 1965–1978 and 1995–2008. The lower inset shows decadal changes in mean annual flow for 10 rivers with outlets in southeast Hudson Bay and James Bay proximal to the study area (highlighted in bold font). Increased flows, particularly in winter, are due, in part, to river regulation during recent decades. The data are from Déry et al. (2011, 2016)

Near the Belcher Islands, located at the southeast corner of Hudson Bay immediately north of James Bay (Fig. 1), the potential for change in oceanography associated with the freshwater balance appears especially acute. In addition to the La Grande River Complex in adjoining James Bay, the duration of the ice-free season has increased by more than a month during the period of 1971 to 2011 in this area (Gagnon and Gough 2005b; Kowal et al. 2017). Furthermore, the region lies within a counterclockwise circulating coastal current (Ingram and Prinsenberg 1998; Saucier et al. 2004), which transports ice melt and runoff collected upstream from the west and south coasts of Hudson Bay. In recent winters, Inuit who live in the community of Sanikiluaq on the Belcher Islands have observed changes in local ice conditions and behavior of leads and polynyas, prompting interest in better understanding the seasonal freshwater balance, the effects of climate change, and the relative roles of sea ice melt and river discharge.

To investigate winter oceanographic conditions and potential interactions between runoff and sea ice melt or brine in southeast Hudson Bay, we initiated a winter study of the coastal waters surrounding the Belcher Islands with a focus on determining the degree of stratification and the relative importance of runoff and sea ice in controlling stratification. Historically, there have been few oceanographic observations in this region and even fewer during winter. The characteristics of the large (> 70 km) La Grande River plume in northeast James Bay have been studied in relation to hydro-related discharge variations in the late 1970s–early 1980s (Freeman et al. 1982; Ingram and Larouche 1987a) and subsequent changes in the plume have been assessed by Hydro Québec (Messier 2002). The characteristics of the much smaller Great Whale River plume under the sea ice along mainland Quebec in southeast Hudson Bay have been studied by Ingram and colleagues (Ingram 1981; Ingram and Larouche 1987b; Ingram et al. 1996). One other study collected profiles of temperature and salinity in waters adjacent to the Belcher Islands in winter (Prinsenberg 1977). Although these studies provided information on winter stratification, none of them collected data in a way that would allow discrimination between runoff and sea ice melt. More recently, Bay-wide oceanographic studies in summer and fall have included several stations located within the southeast corner of Hudson Bay. Of particular relevance to our study, the oceanographic sampling during these campaigns included oxygen stable-isotope measurements for water column profiles, which permit quantitative estimates of the contributions of runoff and sea ice melt (Granskog et al. 2011; Granskog et al. 2007).

Here, we present conductivity-temperature-depth (CTD) profiles, CT time series collected using under-ice moorings, and composition data (salinity and δ18O) for water samples and ice cores collected from the region surrounding the Belcher Islands (Figs. 1 and 2) during four campaigns between late January, 2014, and mid-March, 2015. Tandem measurements of salinity and δ18O were made for the water and ice samples to permit the discrimination between freshwater sources (runoff, sea ice melt). These data are used here to (1) determine the freshwater content of the water column by source (river water, sea ice melt, saline ocean water); (2) establish freshwater inventories for the water column in fall, 2014, early winter, 2015, and late winter, 2015; (3) describe spatial variations in freshwater distribution around the Belcher Islands from early to late winter; and (4) infer temporal changes in the surface waters during the study period using the results of the field campaigns and the composition of sea ice. Based on these results and the past work in this area, and with reference to a NEMO ocean model simulation of winter circulation in the study area, we provide a conceptual model for how river water and sea ice together presently affect the winter oceanographic setting for waters that surround the Belcher Islands.

Fig. 2
figure 2

The icescape in northeast James Bay and southeast Hudson Bay with a fringe of landfast ice along the northeast James Bay coast and surrounding the Belcher Islands and mobile pack ice to the west. A flaw lead, visible as a NNW/SSE trending dark line at the outer edge of the landfast ice, opens intermittently depending on winds. Station locations on the landfast ice platform include Cape Jones (CJ) and 30 surrounding the Belcher Islands as shown in b. Images are modified from NASA’s EOSDIS Worldview on (a) 14 February 2015 and (b) 8 February 2015

Physical Setting

Hudson Bay System

Hudson Bay receives seawater from the Arctic Ocean via the Canadian Arctic Archipelago and exchanges water with the North Atlantic Ocean via Hudson Strait. The circulation at the north end of Hudson Bay between Foxe Basin and Hudson Strait is not well known; however, the primary saline source waters to the Bay are of Pacific origin via Bering Strait (Burt et al. 2016; Jones and Anderson 1994). Within the Bay, the annual mean circulation is robustly cyclonic, transporting much of the freshwater runoff counter-clockwise through the system and then exporting it to southern Hudson Strait (Fig. 1) (Saucier et al. 2004; St-Laurent et al. 2011; Wang et al. 1994). Recent results from a high-resolution ocean general circulation model show that increased river discharge during the spring freshet induces more complex circulation in summer, with multiple small cyclonic and anticyclonic features and mean flow directed through the center of the Bay (Ridenour et al. 2019a). Previous studies identify a coastal or boundary domain, which holds a narrow, swift, river-water rich flow that follows the shore, and an interior Hudson Bay domain, which has slower transport velocities (Saucier et al. 2004; St-Laurent et al. 2011, 2012). The exchange of freshwater between the interior and boundary regions is mainly driven by Ekman transport (St-Laurent et al. 2011, 2012; Ridenour et al. 2019b). In terms of vertical structure, during the open-water season a summer surface mixed layer (SSML) typically occupies the top 30–60 m of the water column, which contains the seasonal freshwater inputs. The surface layer is underlain by cold water extending down to as deep as 125 m in the water column, indicative of the previous winter’s surface mixed layer (WSML) (Granskog et al. 2011). Hudson Bay bottom waters are cold and saline and derived at least in part from overflow of brine-rich bottom waters from Foxe Basin (Defossez et al. 2010).

Although Hudson Bay lies below the Arctic Circle and is ice free during late summer and fall, it becomes completely ice covered in winter (e.g., Hochheim et al. 2011). Its “typical” sea ice cycle involves sea ice forming between September and December (generally from north to south) and then melting between May and August (generally south to north), leaving it ice-free from July/August to November/December (Andrews et al. 2018). This typical ice cycle is changing in response to global warming. Over the period 1980 to 2010, spring sea ice extent decreased 1.4% per year with a ~ 0.3 °C increase in annual air temperature (Hochheim et al. 2011). The total duration of the open-water season increased by 3.1 weeks between 1996–2010 and 1980–1995, with breakup 1.5 weeks earlier and freeze-up 1.6 weeks later (Hochheim and Barber 2014).

Like the Arctic Ocean, Hudson Bay receives large quantities of runoff (> 760 km3 year−1) from an enormous drainage basin (~ 3.7 × 106 km2) (Déry et al. 2011) leading to an annual freshwater yield from runoff for Hudson Bay of approximately 1.0 m (Prinsenberg 1980). An estimated 742 ± 10 km3 of freshwater is stored in sea ice within the bay in April (Landy et al. 2017). To this a further 18 km3 (0.02 m yield) of meteoric water is added annually to Hudson Bay as net precipitation (Ridenour et al. 2019b). With seasonal snow cover typically lasting 4–6 months in central Canada and 6–8 months in Arctic areas, most of the unregulated rivers of Hudson Bay exhibit a nival regime, with low flows in winter when water is stored in the seasonal snowpack, then high flows when the snow melts in spring and early summer. However, winter (January–March) discharge for rivers in eastern Hudson Bay and James Bay exhibits a positive trend over the period 1964–2013, which has been attributed in part to flow regulation for hydropower production during winter months (Déry et al. 2016).

The input of runoff in spring/summer is augmented by meltwater from sea ice, which contributes a freshwater yield of approximately 0.3–0.7 m in northwestern Hudson Bay (Landy et al. 2017), and likely larger quantities in southeastern Hudson Bay (> 1.0 m) as a consequence of the general drift of sea ice southwards during the melt season (Prinsenberg 1988). Freshwater balance within Hudson Bay is ultimately maintained by export of meteoric water and sea ice products (brine or sea ice melt) into Hudson Strait and thence into the Labrador Sea (e.g., Sutcliffe et al. 1983). Although more sea ice is formed in the northwest of Hudson Bay, the tendency for sea ice and sea ice melt to drift southeastward through spring and summer (Landy et al. 2017) means that, annually, sea ice melt exceeds the sea ice produced locally in southeast Hudson Bay (see also Granskog et al. 2011). In winter, the formation of sea ice reverses the melt process by withdrawing similar quantities of freshwater from the surface, leaving behind cold, dense brine.

The Study Area near the Belcher Islands

The Belcher Islands lie about 100 km to the west of Hudson Bay’s east coast and ~ 100 km north of Cape Jones, which is at the mouth of James Bay (Fig. 1). The islands comprise finger-like rock outcrops extending over an area of roughly 3000 km2. The shelf around the Belcher Islands covers approximately 10,000 km2, and water depths vary mostly between 0 and 40 m but exceed 80 m in places. The coastal waters around the Belchers are poorly charted. Freeze-up of waters around the Belcher Islands usually starts in December when air temperatures average − 5 to − 25 °C (Fig. 3), but it typically takes another 6 to 8 weeks for the study area to become fully ice covered. Landfast sea ice accumulates first along the southwestern coasts of the Belchers and later along the north-northeast coasts. By February, landfast sea ice is typically present around the islands, but it is quite narrow, particularly to the west, where there is a large recurrent, north-south trending flaw lead (see arrow in left side of Fig. 2a). To the east of the islands, the sea ice typically remains mobile during January before partial or complete consolidation in February or March. An ice bridge occasionally extends across connecting the Belcher Islands to the eastern coast of Hudson Bay in the vicinity of the Great Whale River (Fig. 1) (CARC 1997). Regional ice analysis charts (Environment Canada) show that average landfast ice thickness is 30–70 cm by February and increases to about 1 m by March, which agrees with our own observations from winter 2015 (Fig. 3b). Break-up usually occurs in late June; between 1971 and 2009, the mean breakup date was 28 June ± 11 days (sd) (Galbraith and Larouche 2011).

Fig. 3
figure 3

Record of temperature (°C) and cumulative freezing °C–day, starting in 1 December 2014 and ending 31 March 2015, at Sanikiluaq (NU). First field campaign in 2015 occurred between January 14 and 31 and second campaign occurred between 17 February and 14 March (highlighted in gray). b Observed ice thicknesses (cm) at 10 stations visited during both of the field campaigns. The line in b shows a second-order polynomial fit (R2 = 0.94, p < 0.001) to the ice thickness vs. age data, consistent with expectations for thermodynamic ice growth (Anderson 1961; Stefan 1889). Points at zero ice thickness reflect dates of freeze-up at these stations as identified by reviewing EOSDIS Worldview satellite imagery

There is no significant streamflow on the Belcher Islands themselves. The nearest river is the Great Whale River located on the mainland, ~ 100 km to the east of the Belchers (Fig. 1). The Great Whale River has mean annual discharge of 19.6 km3 year−1 (1964–2013) with lowest discharge in winter and a pronounced peak during spring snow melt (Hudon 1994). James Bay rivers have a combined mean annual discharge of 270 km3 year−1, of which ~ 110 km3 year−1 is associated with the La Grande complex since completion of the partial diversions of the Eastmain and Opinaca rivers (starting in 1980), upper Caniapiscau River (1982), and Rupert River (2009) (Déry et al. 2016).

Around and among the Belcher Islands, there are numerous small latent-heat polynyas (Smith et al. 1990), which are kept open throughout winter at locations with strong tidal currents (Gilchrist and Robertson 2000). These polynyas and the associated system of flaw leads in southeast Hudson Bay are biologically important (Stirling 1997), providing crucial winter habitat for eider duck populations, polar bears, and seals (Gilchrist and Robertson 2000). Because the Inuit of the Belcher Islands (Sanikiluaq) have relied for generations upon these sea ice habitats including polynyas as sites for harvesting country food (e.g., seals, waterfowl), they have developed a deep knowledge of the seasonal patterns exhibited by the biota that populate these areas in winter. During the 1990s, flaw leads and polynyas around the Belcher Islands began to experience rapid freezing, closing off the over-wintering habitat of the eider duck populations and leading the Inuit to question whether or not these changes were associated with coincident hydroelectric development of rivers entering nearby James Bay. Inuit noted that the seals they hunted started to sink instead of float in winter, which they considered a possible indicator of declines in sea surface salinity or a change in seal diet. Inuit from the community of Sanikiluaq and the Arctic Eider Society approached university researchers to partner on an oceanographic study and improve understanding of the region’s wintertime oceanography.

Material and Methods

Field Sampling

The study was designed and conducted in partnership with the Sanikiluaq Hunters and Trappers Association, facilitated by the Arctic Eider Society. Four field campaigns for sampling were made to the Belcher Islands: (1) 15 January to 5 February, 2014; (2) 27 October 2014; (3) 12–29 January 2015; and (4) 17 February to 14 March 2015 (Table 1). Sampling stations, constrained by the safety of travelling over sea ice, were distributed on the landfast ice platform around the Belchers (Fig. 2). In the winter, the sampling stations were visited by snowmobile. At each station, snow depth and ice thickness was recorded (see Table 2) and a CTD cast was conducted through a hole drilled in the ice. Water samples and ice cores were collected at selected stations (see details in Table 1). In October 2014 (open-water period), four stations were visited by small boat and CTD casts performed and water samples collected. Instrument failure left only bottle samples available from this campaign.

Table 1 Summary of stations and types of samples collected during the four field campaigns starting January, 2014 and ending March, 2015 (W = water sample, CTD = CTD profile, IC = ice core, CT = CT sensor)
Table 2 Freeze-up dates estimated from satellite imagery (MODIS Worldview), observed snow depths and ice thicknesses, and predicted ice thicknesses for January 2015, assuming thermodynamic ice growth as constrained by freeze-up date, observed ice thicknesses, and date of core collection

Eight stations in total were sampled during the first field trip in January 2014 (SK 1, 2, 3, 4, 8, 9, 10, and 14; Fig. 2). Four stations (SK 1, 3, 10, and 9) were sampled in October 2014 (Table 1). In January 2015, CTD profiles were conducted at 17 stations around the Belchers and one station at the mouth of James Bay (CJ; Fig. 2a) and water samples collected at 12 Belcher stations. In February–March 2015, water samples were collected at the same 12 Belcher stations and new stations (T1A, T1B, T1C, A1-A4, and B1-B4) were established in an area southeast of the Belcher Islands and extending 5–20 km towards the mainland (Fig. 2). The ice platform in this area had not been sufficiently consolidated to allow snowmobile travel during the previous campaigns. Ice cores were obtained at a total of seven stations in January–February 2014 and 15 stations in February–March 2015. In January 2015, conductivity and temperature (CT) sensors were moored in the surface waters (depths of 2–3 m as measured from the surface of the landfast ice) at seven stations around the Belchers (SK1, SK3, SK4, SK8, SK14, SK21, SK23) and also at CJ (Fig. 2a). The sensors were left tethered to the ice for periods of 1 to 8 weeks (14 weeks for CJ).

CTD casts were conducted using an RBR sensor during the 2014 field season (temperature accuracy ± 0.002 °C, conductivity accuracy of 0.001 mS cm−1 and pressure (depth) accuracy of ± 0.05%) and a Castaway CTD sensor during the 2015 field season (temperature accuracy ± 0.05 °C, conductivity accuracy of 0.005 mS cm−1 and pressure (depth) accuracy of ± 0.25%). The CT sensors were JFE ALEC Compact-CTs, which have manufacturer’s stated temperature and conductivity accuracy of ± 0.02 °C and ± 0.02 mS cm−1, respectively, and BBRduos, which have manufacturer’s stated temperature and conductivity accuracy of ± 0.002 °C and ± 0.003 mS cm−1.

During winter use of the CTDs, great care was taken to minimize instrument exposure to cold air before the cast. The CTDs were lowered into the water for a two-minute “soak” before beginning a profile, which allows the instrument to thermally equilibrate to the water and the escape of any air bubbles trapped in the conductivity cell. The CTD was then raised to near the surface and then lowered for the profile at consistent speed (about 1 m s−1).

To improve data quality (cf., Halverson et al. 2017), CTD data were post-processed individually. Pressure was corrected for ambient atmospheric pressure as measured at the start and end of the cast while the system was in the air. Any samples where the system was stationary were discarded, and temperature data were reviewed to remove spurious data such as temperature artifacts near the surface (e.g., due to cold instruments). Conductivity data were reviewed to remove erratic measurements near the water surface that can be caused by air bubbles or ice slush trapped in the conductivity flow cell or measurements made when the system is only partially submerged. Historical CTD data (downcasts) collected by S. Prinsenberg, Canadian Hydrographic Service/Canadian Centre for Inland Waters (Burlington), by helicopter (18HE77700) in 1977 were obtained from the Marine Environmental Data Service (MEDS) and are used as-is.

Water samples were collected with a 1.0 L Kemmerer water sampler, which was deployed to designated depths in the water column (measured from the ice bottom in winter) and then brought to the surface to fill dry 1.0 L Nalgene bottles. Sampling depths were chosen to focus on the top 20 m of the water column (1, 5, 10, 20 m) and ~ 5 m above the seabed. In winter, the water samples were stored in an insulated cooler that contained hot water bottles to prevent the samples from freezing during transport. In a field lab established in the Sanikiluaq Hunters and Trappers Association office, the bulk water samples were redistributed into 120-mL Boston rounds and 20-mL borosilicate vials (rinsing three times before filling) for salinity and oxygen isotope analysis respectively, then capped and sealed with parafilm to prevent evaporation and stored at 4 °C.

Ice cores (9.0 cm diameter, maximum of 1 m in length) were collected with a Mark II Kovacs ice coring system. Each core was sectioned immediately by sawing into 5 cm lengths with care to avoid contamination by snow. Individual core sections were placed in sealable plastic bags, stored in a cooler and kept frozen until further processing. At the field lab or University of Manitoba CEOS laboratories, each ice core section was placed in a new, vacuum-sealed plastic bag then left to thaw at room temperature (21 °C) for 24 h. The melted ice samples were then distributed into bottles and vials and stored in the same manner as the water samples.

Lab Analyses

Salinity

The salinities of seawater and melted ice samples were measured with a Guideline Portasal (model 8410A) salinometer in the laboratory of Dr. C. Michel at the Freshwater Institute, Department of Fisheries and Oceans in Winnipeg, Manitoba. The salinometer was calibrated using IAPSO standard seawater of nominal salinity 34.8 before and after each day’s analysis (about 30 samples). A standard deviation of 0.0002 for at least two conductivity readings from a sample was considered to reflect accurate measurement. Based on analyses between June and July 2015, overall accuracy was estimated at ± 0.001 and precision of ± 0.02 (n = 42 duplicates).

δ18O of Water Samples and Melted Ice

Details of the laboratory methods are provided in the SI. Briefly, the oxygen isotopic composition (δ18O) of water and melted ice samples was measured using a Picarro instrument at the University of Manitoba (Walker et al. 2016). Results are reported in per mil relative to Vienna Standard Mean Ocean water (V-SMOW; ‰). Over the course of several months, between June and August 2015, the overall precision obtained for duplicate samples and standards was ± 0.10‰ (n = 53). This precision is similar to that reported in previous studies of coastal Hudson Bay waters (Granskog et al. 2011; Kuzyk et al. 2008).

Analysis of Freshwater Sources

For high-latitude ocean waters where runoff and sea ice interact, tandem measurements of the 18O/16O ratio (expressed as δ18O) and salinity of the water have been developed as a reliable way to determine the relative contribution of river runoff, sea ice melt, and seawater (Ostlund and Hut 1984). More recently, measurements of δ18O and salinity in sea ice cores and water have found application in settings where landfast ice accumulates over ocean water and runoff contribution varies due to advection or plume spreading beneath the landfast ice including the Mackenzie Estuary (Macdonald et al. 1995), Husky Lakes (Macdonald et al. 1999), Lena Estuary (Eicken et al. 2005), Churchill Estuary (Kuzyk et al. 2008), Horton Estuary (Miller et al. 2011), and Svalbard Fjords (Alkire et al. 2015).

Following this previous work, the surface water around the Belcher Islands is considered as saline ocean water (SW) that has been freshened by the addition of runoff, predominantly from rivers (RW), and by the addition of sea ice melt (SIM) during summer. In winter, RW may continue to enter the surface ocean around the Belcher Islands but the formation of sea ice essentially reverses the melting process and withdraws sea ice melt from the ocean, releasing brine. To understand freshwater balance in systems that form sea ice, Östlund and Hut (1984) developed a simple algebraic system to calculate the relative contributions of SW, RW, and SIM to any sample based on the measurement of two conservative water properties (S, δ18O) using the following conservation equations:

$$ {F}_{sw}+{F}_{rw}+{F}_{sim}=1 $$
(1)
$$ {F}_{sw}\ {\mathrm{S}}_{sw}+{F}_{rw}\ {S}_{rw}+{F}_{sim}\ {S}_{sim}={S}_{obs} $$
(2)
$$ {F}_{sw}\ {\updelta}^{18}{\mathrm{O}}_{sw}+{F}_{rw}\ {\updelta}^{18}{\mathrm{O}}_{rw}+{F}_{sim}\ {\updelta}^{18}{\mathrm{O}}_{sim}={\updelta}^{18}{\mathrm{O}}_{obs} $$
(3)

where F refers to the fraction of each water type in the given sample denoted by the subscripts sw, rw, and sim. S and δ18O refer to the salinity and isotopic composition of the three components making up the water (subscripted sw, rw, sim) or the observed values of the given water sample (subscripted obs). For the Belcher Island study, this three-component system (Fsw, Frw, Fsim) was solved quantitatively for each seawater sample after assigning values appropriate for the local setting for S and δ18O for each of the three water types (or end-members).

Interpreting δ18O Records in the Sea Ice

The chemical properties of sea ice indirectly reflect the properties of the surface waters from which the ice was grown. In the case of salt (salinity) and other dissolved components, most of which are rejected as brine during freezing, the sea ice contains a poor, non-conservative record of the seawater. In contrast, the δ18O composition of the sea ice reflects that of surface ocean water but with a fractionation offset such that the sea ice is isotopically heavier than the water it formed from. Thus, the δ18O profiles in the ice cores collected in this study may be used to infer changes in the meteoric water content of surface waters throughout the period of ice growth. To make this inference, we require (1) a quantitative estimate of the fractionation during freezing and (2) a relationship between depth within the ice core and time of freezing.

Estimate of Fractionation during Freezing

Fractionation for formation of Belcher Island sea ice (2.1 ± 0.1‰) was determined empirically based on tandem measurements (n = 20) of the δ18O in ice near the bottom of the core and seawater below the ice at time of core collection (cf., Kuzyk et al. 2008). This measured fractionation falls within the range of 2–3‰ established in other studies (Crabeck et al. 2014; Eicken 1998; Macdonald et al. 2002).

Estimate of the Depth-Time Relationships Within the Ice Cores

The conversion of depth in sea ice cores to time of freezing was based partly on the methods used by of Macdonald et al. (1995). First, each ice core was constrained by applying two known dates: (i) initiation of freeze-up (top of the core) estimated from MODIS satellite imagery (NASA, EODIS Worldview) and (ii) the date of core collection (bottom of the core) (Table 2). These two constraining dates, together with any additional dates when sea ice thicknesses were measured (see Fig. 3b), were used to relate depth in the core to time using Stefan’s equation for ice growth (Anderson 1961; Stefan 1889). Stefan’s equation provides an analytical solution for thermodynamic ice growth as a function of freezing-degree days (t) in the form of:

$$ H=d\ {t}^{0.5}, $$
(4)

where H is the ice thickness at time, t, and d is a constant that varies among cores. Comparison of dates developed for ten cores with ice-thicknesses measured 35–53 days into the period of ice growth (Table 2) suggested that dates assigned to sections in the middle of the cores, where errors would be greatest, could be out by as much as 1 week. Detailed explanation of the procedure for converting depth in ice cores to time of freezing is provided in the SI.

Simulation of Circulation in Southeast Hudson Bay in Winter

A 3-D numerical ocean model (NEMO) based on the Nucleus for European Modelling of the Ocean version 3.4 and coupled to the Louvain-la-neuve Ice Model version 2 was used to simulate circulation in Hudson Bay, including the periods of intensive data collection (January–March 2014 and January–March 2015). Within the BaySys project, a Bay-wide initiative to investigate effects of hydroelectric regulation and climate change, this model has recently been used to examine the sensitivity of Hudson Bay’s freshwater dynamics to runoff forcing and to contrast spring/summer vs. fall circulation patterns (Ridenour et al. 2019a, b). For the simulation in this study, the model configuration was similar to that in the works cited above. Specifically, we used the Arctic and Northern Hemisphere Atlantic with 1/12° resolution (ANHA12) (Hu et al. 2018), which yields a horizontal resolution of 3.5–5.5 km for the study area. The vertical has 50 geopotential levels with the highest resolution (∼ 1 m) in the top 10 m. The initial fields for salinity, temperature, horizontal velocities, etc., as well as boundary conditions, were from GLobal Ocean ReanalYsis and Simulations (GLORYS2v3), and the atmospheric data to force the model was from the Canadian Meteorological Centre. The river discharge in Hudson Bay is monthly and is based on the Dai et al. (2009) dataset. For the purposes of this study, we show the modeled sea surface height and calculated surface geostrophic velocities averaged over a 3-month period (January, February, and March for each of 2014 and 2015). For ease of comparison with previous work (Ridenour et al. 2019a), sea-surface heights are referenced to height above the geoid (absolute dynamic topography), which may be positive or negative.

Results

Spatial and Temporal Variations in Seawater Salinity and Temperature

Based on CTD profiles collected around the Belcher Islands, the salinity ranged from about 26 to 30.5 and the temperature from − 1.5 to − 0.5 °C during the 2015 winter study (Fig. 4). Similar salinities and temperatures were observed in the study area during January–February 2014 (Petrusevich et al. 2018). CTD profiles at a shallow site (12 m) near Cape Jones at the mouth of James Bay in January 2015 showed uniform salinity of about 25 and temperature of about − 1.3 °C (not shown but see CT record in Fig. 6).

Fig. 4
figure 4

Salinity and temperature profiles for the shallow waters surrounding the Belchers for January 2015 (a) and February–March 2015 (b) shaded according to latitude. Stations mentioned in the text are highlighted

The common features in the Belcher Islands CTD profiles in both January 2015 (Fig. 4a) and February–March 2015 (Fig. 4b) are (i) a fresher surface layer with lower salinities towards the south/southeast (e.g., SK 1, SK 3) and (ii) the presence of warmer, salty water in the subsurface below 30–40 m water depth (e.g., SK 2). The warm, saline deep water was interpreted in previous studies as a regional feature established prior to freeze up (Petrusevich et al. 2018). In January 2015, the subsurface waters below a depth of 40 m around the Belchers had a salinity of at least 28.4 and temperature varying from − 1.2 °C to as high as − 0.6 °C (Fig. 4a). In February–March 2015, the salinity of the subsurface waters around the Belchers had increased slightly to between 29 and 30.5 and the temperature was between − 1 and − 0.4 °C (Fig. 4b).

The fresher surface layer and associated shallow stratification of the water column was observed at most stations south/southeast of the Belchers during winter 2015 (Fig. 4 and Fig. SI-1). In January 2015, stations SK 1 and SK 3 exhibited the lowest surface salinities of the Belcher Islands study area (25.8–26). In February–March 2015, the cluster of stations T1A, T1B, T1C, A1-A4, and B1-B4, located 5–20 km off the southeast coast (Fig. 4), had the lowest surface salinities (< 26.5) (Fig. 4b; Fig. SI-1). These stations were sampled for the first time in February–March 2015 due to delayed ice consolidation in this area in 2015. The salinity profiles for the most southeasterly stations (T1A, T1B, T1C, A1-A4, and B1-B4) show a stratified water column with fresher surface water extending to approximately 20 m in depth (Fig. 4b). A few southerly stations such as SK2, located among the islands rather than around the perimeter (Fig. 4), do not show as fresh a surface layer as stations SK1 and SK3, despite similar latitude (Fig. 4a, b). This difference may be caused by restricted flow through the channels. The strongest contrasts to southern stations, however, are northern and western stations such as SK4 and SK14 (highlighted on Fig. 4a, b), which show saline surface waters and no stratification. These cold, saline, well-mixed profiles match expectations for winter surface mixed layers formed by brine-driven winter convection in the absence of winter river inflow (Granskog et al. 2011; Saucier et al. 2004).

The spatial distribution of salinity and temperature during each sampling period is illustrated as a sectional path around the Belcher Islands beginning at the northwestern station SK 12, extending counter-clockwise around the islands, and terminating back at the northeast corner at station SK 15 (Fig. 5). The stratified water column and fresher, near freezing, surface waters occur to the south/southeast of the Belcher Islands in both January (Fig. 5a, b) and February–March 2015 (Fig. 5c, d). In January 2015, a lobe of low salinity surface water was evident along the southeast coast of the islands (SK 8 to SK 23) with the lowest surface salinities (~26) at SK 1 (Fig. 5a). In February–March 2015, the low salinity surface layer appears to have expanded laterally (westward and northward) and deepened (Fig. 5c) with the layer attaining depths > 20 m at some stations (see e.g., SK 1 in Fig. 5c). Stations SK8 and SK23 are more affected by the fresher surface waters in February–March 2015 compared to January 2015. The deep waters in February–March 2015 have a salinity of 29–30 and a temperature near 0.5 °C, which is well above freezing (− 1.6 °C). A second section orthogonal to the southeast coast of the Belcher Islands that incorporates stations SK22, T1C, T1B, and T1A (see arrow on lower right on the inset for Fig. 5) shows a strongly stratified water column over the 20 km distance during February–March 2015 with surface salinity ~ 26.5 (Fig. 5e).

Fig. 5
figure 5

Sectional view of salinity and temperature for January 2015 (a, b) and February–March 2015 (c, d and e, f), using in a, b and c, d, a section that extends around the Belcher Islands beginning in the northwest (station SK 12), extending down the west side of the islands, then across the south side, and finally up the east, and in e, f, an orthogonal section extending from the southeast coast of the Belcher Islands eastward ≈ 20 km towards mainland Hudson Bay (see arrow on inset map)

The continuous CT records from ice-tethered moorings placed at 1 m distance beneath the sea-ice bottom (depth of ~2–3 m) show time variations in salinity between January and February–March 2015 unique to particular stations or subsets of stations (Fig. 6). The CT sensor at Cape Jones (CJ; Fig. 6) recorded an average (range) salinity of 24.7 (22.8–25.8) over the period 26 January 2015–31 March 2015 with no particular temporal trend (Fig. 6). Around the Belchers, sensors placed at south/southeast stations such as SK1 and SK3 show low salinities in January, which increase into February–March 2015. At SK1, a salinity of 25.5 was recorded on 16 January 2015; salinity increased abruptly from 16 to 25 January, and thereafter salinity averaged about 26.7 (Fig. 6). At SK3, the CT record shows an average surface salinity of 26.2 during the period 16–28 January, followed by an almost linear increase in salinity until late February, at which time the surface salinity was 27.1 (Fig. 6). In contrast to SK1 and SK3, stations SK8 and SK4 show initially high salinity of 28.7–29 followed by an abrupt decrease to about 28 on 25 January 2015 (Fig. 6). The decreasing salinity at these two locations is consistent with the low salinity surface layer having expanded horizontally and spread westward between January 2015 and February–March 2015. SK23 shows relatively stable salinity (~ 27.5), slightly fresher than station SK8. SK21 is slightly fresher again with a pronounced semi-diurnal tidal signal. At the most westerly station (SK14), salinity is about 29.4 in January 2015 and gradually increases throughout February to a maximum of about 29.75 in early March.

Fig. 6
figure 6

Salinity time series collected using under-ice moorings during winter 2015. See Fig. 2 for station locations. CJ refers to the station at Cape Jones at the mouth of James Bay and the other stations surround the Belcher Islands

A comparison of water properties at stations SK 1 and SK 3 in January 2015 (Fig. 4a) with those during the preceding October 2014 (data shown in Table SI.1, Supplemental Information) reveals that surface salinity decreased between fall and winter (from 27.7 at SK 1 and 28.0 at SK 3 to ~ 26 at both stations), while deep water salinity increased (from 28.9 to 29.8). Unfortunately, data from October 2014 are not available for stations on the north or west side of Belcher Islands. However, in October 2005, the surface salinity west of the Belcher Islands (station BI-2) averaged 27.5 (ArcticNet station AN02; unpublished data not shown). Thus, it appears that surface waters are fresher in winter than fall to the south/southeast of the Belcher Islands, but more saline in winter compared to the fall to the west of the Belcher Islands.

Oxygen Isotopic Ratios (δ18O) in Seawater and River Water

The seawater samples collected from the waters surrounding the Belcher Islands during 2014–2015 exhibited δ18O values between − 5.2 and − 2.3‰ (Fig. 7; Table SI.1). These values lie within the range previously established during the open water season for southeast Hudson Bay (− 5.7 to − 2.3‰; ArcticNet/CCGS Amundsen cruises). Rivers discharging into western Hudson Bay and James Bay “upstream” from the Belcher Islands study area are reported to have δ18O values varying from − 19.5 to − 10.3‰ (Fig. 7 data along the y-axis (Granskog et al. 2011)). River discharge from La Grande Complex sampled in January 2015 had a δ18O value of − 13.28‰ (n = 1, unpublished data).

Fig. 7
figure 7

Plots of δ18O versus salinity for the study region. a Hudson Bay water (study region water samples (plus signs) and ArcticNet/Amundsen (AN) cruises 2005–2010 (empty squares)), river samples (empty diamonds), and melted sea-ice samples (empty circles). End-members assigned for local seawater, river runoff, and sea ice (Table 3) are shown as (filled diamonds). The line indicates the mixing line between river runoff and seawater end-member values. b Seawater samples from this study replotted to highlight the changes in properties between January–February (plus sign), October 2014 (cross), January 2015 (filled circles) and February–March 2015 (empty square). Values above the mixing line indicate the presence of sea-ice melt (+SIM), while values below the mixing line indicate the presence of brine (–SIM)

When δ18O data are plotted against salinity for the samples collected during our study and previously published data, three groupings are evident (Fig. 7; Table SI.1). Seawater samples, clustered toward the top right in Fig. 7, show a relatively large component of variance along a mixing line between seawater (farthest right) and river water (bottom left). However, seasonal shifts in δ18O-S water properties are also evident among the four campaigns as overlapping clusters of data where the first campaign spanning January–February 2014 sets the lower δ18O limit for out dataset (most negative) and the second campaign in October 2014 generally sets the upper δ18O limit (Fig. 7b). The winter campaigns in 2015 lie between the limits, displaying a slight shift to more negative δ18O values from January 2015 to February–March 2015. When the δ18O data for surface water samples are examined, there is a clear spatial pattern in δ18O values during the winter 2015 campaigns, with lowest δ18O values (below − 4‰) in surface waters southeast of the Belchers and highest values to the north and west of the islands (Fig. SI-2).

Oxygen Isotopic Ratios (δ18O) in Sea Ice

Sea ice samples collected during January–February 2014 had salinity and δ18O values varying from 4 to 14 and − 4.0 to − 1.0‰, respectively. Sea ice samples collected during February–March 2015 had δ18O values similar to those in 2014 (p > 0.05) varying from − 3.5 to − 0.5‰ (salinity not available). The 2014 ice samples plot in the upper left in Fig. 7, reflecting their low salinity relative to seawater samples of similar isotopic values. Across all the sea ice samples from both years, δ18O values increase significantly with increasing latitude (r2 = 0.18, p < 0.001, n = 330; Fig. SI-3). Thus, the spatial pattern in δ18O in sea ice roughly reflects the spatial pattern in δ18O in sea water described above.

Vertical profiles of δ18O in ten sea ice cores collected during February–March 2015 from sites along the transect around the perimeter of the Belcher Islands (see map inset in Fig. 5) are shown in Fig. 8. Within most ice cores, there are one or two samples near the top with anomalously low δ18O values (open circles in Fig. 8). Isotopically light values in the top 10–15 cm of the sea ice, which represents the earliest stage of sea ice formation, likely reflect isotopically light snowfall (about − 25‰ in southern Hudson Bay (Smith et al. 2015)) incorporated into the sea ice as it first began to form. Disregarding these anomalous surface samples, there are statistically significant decreases in δ18O with depth in the ice at SK1, SK21, SK22, and SK23, with the decreases mostly happening step-wise at depths between 20 and 30 cm (Fig. 8), or in the case of SK21, over two steps. Other cores, including SK4, SK10, SK14, and SK15, show statistically significant increases in δ18O with depth, with the increases generally occurring most clearly in the upper (earlier) portions of the profiles (Fig. 8). Note that these upper portions of the cores (top 30 cm) mostly reflect ice growth that occurred before our January 2015 field campaign as ice thicknesses observed at that time varied from 28 to 90 cm (average 50 cm; Table 2).

Fig. 8
figure 8

Individual vertical profiles of δ18O in ice cores collected along the section extending around the Belcher Islands (see map Fig. 5) during the 2015 winter field season. Increasing depth in the ice reflects ice accumulated later in the year such that change in ice composition with depth reflects change over time in the composition of surface water from which the ice was grown. The top of the core (0 cm) represents the start of ice formation at the site; the bottom of the core reflects the last few days up until the ice core was collected. Low δ18O values in the top 10 cm (open symbols) are likely due to contamination by snow

Modeled Circulation

Simulated surface circulation for Hudson Bay for January–March 2014 and 2015 shows that the Belcher Islands study area lies directly in the path of generally northward-flowing James Bay outflow (Fig. 9a,b). To the west of the Belchers, during both years, there is very strong northward flow through the center of Hudson Bay associated with the west Hudson Bay cyclonic circulation cell (Fig. 9). In January–March 2014, the relatively buoyant James Bay outflow produces strong surface flow (> 0.05 m s−1) west of the Belcher Islands driven by east-west sea-surface height gradients (Fig. 9a). Sea-surface heights (referenced to the geoid) are high throughout eastern James Bay and all along the east coast of Hudson Bay. The strong northward-flowing boundary current simulated for January–March 2014 is not dissimilar to the mean fall circulation in Hudson Bay previously shown (Ridenour et al. 2019a). In January–March 2015, the simulation shows James Bay outflow mostly heading northeast when it enters Hudson Bay but then some westward deflection around the south side of the Belchers (Fig. 9b).

Fig. 9
figure 9

NEMO ocean model simulations of winter circulation in the study area during January–March 2014 (a) and January–March 2015 (b). Arrows represent surface geostrophic current velocities and colors reflect sea surface height (absolute dynamic topography or height above the geoid)

Discussion

Distributions of Freshwater in Winter

Lacking winter oceanographic data or annual mooring records for southeast Hudson Bay due to its remoteness and inaccessibility by ship, we expected to observe, at all sites, a continuous increase in the salinity of surface waters reflecting the addition of brine from the growing ice cover and the deepening of the pycnocline throughout winter because this is the seasonal trend observed in west Hudson Bay (Prinsenberg 1987). However, two distinct oceanographic domains are evident in our dataset based on the spatial variability in salinity and temperature during the 2015 (Figs. 4, 5, and 6) and 2014 (Petrusevich et al. 2018) winter campaigns. In the east-southeast part of the study area, the water column was strongly stratified throughout winter with an approximately 20 m thick, relatively fresh surface layer. The west-northwest part of the study area had a water column that was well-mixed and more saline and became increasingly saline as winter progressed, which is generally similar to the seasonal variations in surface properties observed in west Hudson Bay. The two domains are present in the winters of both 2014 and 2015, but we focus our discussion below on 2015 because we obtained a larger dataset in that year, augmented with continuous records of the salinity just beneath the ice (Fig. 6).

The lobe of freshwater constrained by the southeast coast of the Belcher Islands is sharply defined in the January section but appears to have expanded in area and depth by late February (compare Fig. 5a and c). The 20 km section out from the coast, collected in February 2015 (Fig. 5e, f), shows that the stratification and relatively fresh surface layer at the southeast corner of the Belcher Islands extended well into the strait between the island and the mainland to the east. Beneath this freshwater lens south/southeast of the Belchers, the seawater was noticeably warmer and well above its freezing point, especially in late February–March 2015 (Fig. 4). This warm seawater is presumably a remnant from summer that has remained isolated from heat exchange with the cold atmosphere due to stratification (Petrusevich et al. 2018; Prinsenberg 1984) and is evidence that typical winter convection has not occurred in this area. Waters in many regions of the Arctic remain stratified year-round (cf., Peralta-Ferriz and Woodgate 2015), and if the halocline is strong enough, trap summer heat in the shallow subsurface throughout the winter (Jackson et al. 2010). To the north and west of the Belcher Islands, where stratification was weaker or absent, the surface waters cooled and became slightly saltier between January and February/March (Figs. 4, 5, and 6), consistent with the seasonal trend that has been observed in west Hudson Bay (Prinsenberg 1987).

Broadly speaking, the separation between fresher, well-stratified seawater to the east/southeast of the Belchers and more saline, less-stratified seawater to the north/west of the Belcher Islands makes sense because the islands lie at the outer boundary of the coastal current in southeast Hudson Bay, which transports much of the runoff around the margin of Hudson Bay and eventually out into Hudson Strait (e.g., St-Laurent et al. 2011). Some authors have assumed that the relatively fresh coastal current behaves like a conduit, while others have shown that perhaps 25% of the runoff is diverted into the interior of Hudson Bay via Ekman transport in summer, and then released back to the boundary area in fall (St-Laurent et al. 2011, 2012). Regardless, the simulated circulation (Fig. 9), together with our observations (Figs. 4, 5, and 6), indicate for the first time that the distinction between a fresher nearshore region and a saltier interior ocean, previously observed during summer (cf., Granskog et al. 2007; Lapoussière et al. 2013), is maintained throughout winter. The roles of sea ice and runoff in maintaining this structure remain unclear until we apply salinity-δ18O tracer pairs to distinguish the runoff from processes associated with the sea ice (melt/formation).

Application of S and δ18O to Determine Contributions of Fresh Water Sources

The spread in salinity and δ18O measured in the water samples (Fig. 7) shows that both SIM and RW are important additions to seawater (SW) in the water column. River water, which has been distilled through atmospheric processes, exhibits low δ18O values ranging from − 10 to − 20‰ (Fig. 7a, bottom left). The seawater samples collected around the Belcher Islands during our study, and prior to it, form a group in the S range 26 to 30, and δ18O range − 5 to − 2‰ (Fig. 7a, top right). These data cluster around a major axis of variance represented by the black line in Fig. 7, which is produced by mixing saline, isotopically-heavy ocean water with isotopically light river water. As found in numerous other studies where sea ice is a seasonal feature, there is a second component of variance orthogonal to this line. This second source of variance is produced by the melting and formation of sea ice (+SIM and −SIM, respectively). In winter, sea ice formation expels brine into the water with the effect (in Fig. 7) of shifting data points towards higher salinity and slightly more negative δ18O. In summer, melting of sea ice adds water of low salinity and relatively high δ18O as shown by the cluster of points at the top left in Fig. 7a.

The algebraic approach to determining relative contributions of RW, SIM, and SW to seawater samples, pioneered by Ostlund and Hut (1984), has been successfully applied at many locations in the Arctic Ocean (Alkire and Trefry 2006; Bauch et al. 2011; Guay et al. 2001; Macdonald et al. 1995; Newton et al. 2013; Yamamoto-Kawai et al. 2008) and in Hudson Bay (Granskog et al. 2009, 2011; Kuzyk et al. 2008; Tan and Strain 1996). The algebraic method has limitations: two conservative tracers (S, δ18O) and the conservation constraint, described in the “Material and Methods” section, permit solutions for only three water-mass components (RW, SIM, SW). Furthermore, adding brine to seawater is equivalent to forming sea ice and thus removing SIM from the water, which means that sea ice formation decreases SIM values and negative values for SIM are equivalent to net brine production. A crucial step to producing algebraic solutions for the water samples is the assignment of appropriate end-member properties (δ18O, S) for the three water-mass components that make up surface water in a given study area (RW, SIM, SW). Study objectives and geographic settings vary widely and, therefore, so do the end-member values chosen for water mass studies. Taking into consideration (1) the objective of this study to investigate freshwater content and sources in surface water (~ 40 m) throughout winter around the Belcher Islands, (2) the dataset collected during this study (water and sea ice), and (3) previous sampling of Hudson Bay water for salinity and δ18O (e.g., Granskog et al. 2009, 2011), we have chosen the end-member values listed in Table 3.

Table 3 Local end-members assigned for this study and mean (SD) fractionation factor as determined from ice and under-ice water sample pairs (n = 20)

For the RW end-member (− 13.3‰ for δ18Orw), we have used the flow-weighted average δ18O for measurements available for rivers upstream (along the Hudson Bay coast to the west and north and in James Bay) of the Belcher Islands (Table 4). We also included the Little Whale River (LWR) and Great Whale River (GWR), which discharge into southeast Hudson Bay (Fig. 1) to the east of the Belchers, because they are proximal to the study area. However, their low annual discharge (3.7 km3 year−1 and 19.8 km3 year−1, respectively), and particularly their low discharge in winter (e.g., 3.2–9.5 km3 year−1 for the GWR (Hudon et al. 1996)), indicate they have little potential to influence water around the Belchers. For comparison, the mean annual discharge of the La Grande regularly exceeds 100 km3 year−1 (Déry et al. 2016) and has peak flow exceeding 126 km3 year−1 between December and February (Messier 2002). The LWR and GWR are also located > 75 km away from the Belchers, across a channel that, in winter, contains mostly pack ice. According to Ingram and Larouche (1987b), who directly measured the extent of the under-ice plume of the GWR in relation to ice cover and river discharge, the pack ice in this area limits plume expansion westward. These workers found that the GWR plume was deflected northward from the river mouth under the landfast ice cover. Because of the presence of landfast ice, which suppresses wind mixing, the plume was larger in winter under the ice cover compared to open water periods with similar discharge. However, according to their Fig. 2, the 25 psu surface isohaline extends at most 30 km away from the coast under the landfast ice (Ingram and Larouche 1987b). The plume would have to be more than double this extent to influence the Belchers. The NEMO simulation for January–March 2015 shows complex, possibly anticyclonic circulation in the channel between the Belcher Islands and the mainland to the east (Fig. 9), which seems to contradict the observations of northward flow of the GWR plume (Ingram and Larouche 1987b) but the NEMO model used in this study does not well resolve small individual river plumes nor circulation very near to the coast. The LWR has only 19% of the discharge of the GWR and is further removed from the Belchers in a northeasterly (downstream) direction, which means it has even less opportunity to supply river runoff to the Belchers in winter. Regardless of the extent of influence of the LWR and GWR, their addition or removal makes a negligible difference to the flow-weighted average (Table 4) used for RW end-member assignment.

Table 4 Mean annual flow and δ18O values for rivers upstream or proximal to the study region, organized by territory, province, and drainage region (see Fig. 1 for locations). The mean annual flow data are from Déry et al. (2011, 2016) and Stadnyk et al. (2019), the La Grande δ18O data are from this study (winter 2015), and the remaining δ18O data are from Granskog et al. (2011) = and/or 2005–2010 ArcticNet/CCGS Amundsen cruises (unpublished)

In contrast to GWR and LWR, the La Grande River δ18O value exerts a relatively strong influence (roughly 30% or 110 km3 year−1/362 km3 year−1; Table 4) on the flow-weighted average δ18O value used for RW end-member assignment. This heavy weighting arises because we have no data for other rivers discharging into James Bay. However, fortuitously, this weighting (30%) reasonably represents the contribution of all James Bay rivers to total Hudson Bay discharge (i.e., 35%: 270 km3 year−1/760 km3 year−1) and the δ18O data for the La Grande (− 13.28‰) is a reasonable estimate of the average δ18O value for these rivers. Granskog et al. (2011) inferred a mean value of − 13.0 ± 0.3‰ for all rivers discharging into James Bay (see their Table 2) based on using the zero intercept for salinity–δ18O data for an ocean section crossing the mouth of James Bay. The value of − 13.0‰ also is consistent with the estimate for precipitation in the James Bay watershed (Gibson et al. 2005). Thus, the RW end-member assigned in this study appears to be regionally representative despite the scarcity of δ18O data for rivers in some subregions.

Net precipitation (P-E) also adds (or subtracts) distilled (isotopically light) freshwater to the sea surface. For the southern Hudson Bay region, the isotopic composition of precipitation has been measured or estimated at − 12 to − 15‰ (Delavau et al. 2011; Gibson et al. 2005). In view of the small (~ 0.02 m) and uncertain contribution of P-E relative to the riverine yield for Hudson Bay (Prinsenberg 1977; Ridenour et al. 2019b; St-Laurent et al. 2011), and its similarity in isotopic composition to the upstream river waters, we henceforth consider RW to include both runoff and precipitation, the bulk of which will be runoff. An uncertainty of ± 1.3‰ (SD) assigned to RW in Table 3 includes variation associated with δ18O in river water and in the mean flows of the rivers. Small and unregulated rivers have an inherently larger uncertainty because of pronounced seasonal variability in δ18O (cf., Kuzyk et al. 2008), but individually they contribute less to the total inflow. δ18O values are available for the two largest upstream river systems, the Nelson and the La Grande Rivers. Both are regulated river systems containing lakes and reservoirs. The Nelson has a relatively constant flow rate over the annual cycle with small standard deviation (~ 6%). Based on a large dataset (n = 344) obtained from upstream sites on the Nelson River, Smith et al. (2015) estimated the coefficient of variation as ± 4% of the average δ18O value (− 10.64‰). The variability in the La Grande River δ18O value is not well known. Arctic rivers generally contain a slight amount of salt (see, e.g., Alkire et al. 2017; Macdonald et al. 1995), which we have also allowed for in Table 3.

The seawater end-members (SW) are based on the S-δ18O data pairs collected in this study together with data collected during the ArcticNet 2005–2010 ship-based sampling at stations nearby to those in the study region (Granskog et al. 2011; unpublished). The salinity and δ18O end-member values of 32.2 ± 0.1 and − 2.7 ± 0.1‰, respectively, lie to the right of our data set as shown in Fig. 7a, and on the apparent mixing line on data collected during ArcticNet cruises.

For the sea-ice δ18O end-member (SIM), we use a value of − 0.7‰, which is consistent with the SW end-member and the estimated average fractionation of 2.1 ± 0.1‰ based on coincident measurements of sea ice and seawater (n = 20). This value falls within the range of 2.0–2.2‰ estimated by Kuzyk et al. (2008) and Granskog et al. (2011) for other locations in Hudson Bay. Salinity is not conservative in sea ice, and tends to vary over time and with depth in ice cores (e.g., Ehn et al. 2007; Leppäranta 1993; Macdonald et al. 1995). Based on the ice cores analyzed in this study, we use the average sea ice salinity of 6.1 ± 1.4 as the SIM salinity end-member.

Varying the δ18O and salinity end-members within their estimated errors (Table 3) changes our estimated proportions of river water by as much as 0.016 and our estimated proportions of sea ice melt by as much as 0.019. Comparisons among stations and times within this study are much better than implied by these errors because the algebraic solutions shift for all samples ensuring good internal consistency; however, apparent differences in RW or SIM fractions of < 0.035 between this and other studies (using different end-members) should be viewed with caution.

Distributions of RW in Winter

The partitioning of the freshwater content in the seawater samples into the RW and SIM components (Figs. 10 and 11) shows that RW dominated the freshwater contained in the 20 m surface layer at the southeast corner of the Belcher Islands in January–March 2015. The fraction of RW (Frw) in the upper 20 m in this area exceeded 0.125 in January 2015 and 0.15 in February–March 2015 (Fig. 10). The RW was present predominantly within a lobe abutting the southeast shore of the Belcher Islands; the lobe widened and deepened between January and late February. The inventories of RW in the top 20 m, shown as bars across the top of the sections (Fig. 10 and see Table 5), were > 2.0 m at stations SK1, SK21, SK3, and SK23 in both January and February–March 2015. As an indication of the uncertainty in the inventory that could be produced by varying the RW end-member within the bounds given in Table 4 (− 12 ± 2‰) (corresponding inventories for the top 20 m in January at SK1 would be 2.2 ± 0.2 m. Similar large (> 2.0 m) inventories of RW were present in the top 20 m of the water column along the transect orthogonal to the southeast coast of the Belcher Islands (see stations SK 22, T1C, T1B, T1A; Table 5). Thus, the waters are similarly rich in RW in February–March 2015 within ~ 20 km of the southeast tip of the Belcher Islands.

Fig. 10
figure 10

Sectional view of the fraction of river water (Frw) present in the water column along the section extending around the Belcher Islands (see map Fig. 5) in January 2015 and February–March 2015. Above the sections, bars show the RW inventory in the top 20 m of the water column. For reference, red bars show the mean (± SD) RW present in the top 20 m in October 2014 and January 2015

Fig. 11
figure 11

Sectional view of the fraction of sea-ice melt (Fsim) present in the water column along the sectional view (see map Fig. 5) of the study area in January 2015 (a) and February–March, 2015 (b). Above the sections, bars show sea ice thicknesses and the SIM inventories in the top 20 m of the water column. For reference, red bars show the mean (± SD) SIM present in the top 20 m in October 2014 and January 2015. Note that negative SIM inventories imply net brine

Table 5 Inventories of river water (RW) and sea ice melt (SIM) in the top 20 m of the water column during the four field campaigns

On the northwestern side of the Belchers at stations SK12, SK4, and SK14, the proportion of RW (Frw) in the upper 20 m was < 0.025 in January and < 0.075 in February–March (Fig. 10). The inventory of RW in the top 20 m was relatively small in January (~ 0.5 m), compared to the > 2.0 m found in the southeast lobe. At SK 4 in the northwest corner, the inventories (Table 5) suggest that > 1 m of RW was added to the top 20 m after January. The increase in RW fraction at SK 4 extended down to the seafloor (Fig. 10) with the resultant total water column RW inventory increasing by over 4 m (Table 6). These inventories underscore the far stronger RW presence in the waters southeast of the Belchers compared to the northwest, but also show that RW was being transported westward (e.g., SK 4) episodically during winter. The late-January decrease in salinity under the ice at SK4 and SK8 recorded by the CT sensors (Fig. 6) may represent one of these episodes. In the δ18O profile in the ice core from station SK4, there are several small decreases in an otherwise increasing downward trend (e.g., at 35 cm, 65 cm, 75 cm; Fig. 8), which also may reflect episodic intrusions of RW under growing sea ice at a site otherwise becoming saltier due to the continued release of brine.

Table 6 Inventories of river water (RW) and sea ice melt (SIM) in the total water column during the four field campaigns

As shown by Macdonald et al. (1995) and Macdonald and Carmack (1991), landfast ice may record variations in surface-water δ18O properties as the ice thickens, but with the ice having an approximately + 2‰ offset in δ18O due to isotopic fractionation during freezing. To exploit this record, the fractionation must be accounted for and, additionally, the depth in the ice must be transformed to time (see details in the “Material and Methods”). Figure 12 shows temporal variations in surface-water Frw for the period between freeze up (mid-December 2014) up to the time the ice cores were collected (mid-March 2015) as derived from the ice-core records using the techniques of Macdonald et al. (1991, 1995). The same two domains seen in the water-column sections—northwest vs. southeast of the Belchers—are clearly evident in the ice core record (Fig. 12), showing that the division into two oceanographic domains is a persistent feature of the Belcher Island coast throughout winter including the intervals prior to and between our sampling campaigns. Northwest of the Belchers, the ice core record at station SK 14 shows Frw ~ 0.06 soon after freeze-up in December, a slight decrease as winter and ice growth progressed to a minimum of about 0.01 in January, and then a return to about 0.06 in late February-early March. The under-ice CT at SK14 also captured the slight freshening between January and late February-early March (Fig. 6) but obviously lacked the information about December. Southeast of the Belchers, the ice core record at SK21 and SK22 shows that the greatest proportions of river water (Frw > 0.15) preceded our January water sampling. Relatively high fractions of river water (Frw ~ 0.10–0.15) remained in the surface water southeast of the Belchers throughout the rest of winter according to the ice core record (Fig. 12), in agreement with the water-mass distributions in late February and March 2015 displayed in Fig. 10.

Fig. 12
figure 12

Record in the sea ice along the section extending around the Belcher Islands (see map Fig. 5), depicted as the fraction of river water (Frw) of the surface waters from which the ice was grown. The vertical axis is time beginning at first freeze-up date (10 December 2014) and ending last day of sample collection (14 March 2015)

The sparser data set for January–February 2014 (Tables 5 and 6 (not plotted)) show generally the same spatial pattern seen in winter, 2015. For the entire water column, which varied in depth among stations (Table 6), the storage of RW in winter 2014 was as high as 9.7 m within the southeast lobe (SK1) compared to 1.9–7.4 m in the northwest domain (SK4, SK8, SK14). Although the spatial patterns are similar between the two winters, there was significantly more RW to the northwest of the Belchers in 2014.

Distributions of SIM in Winter

Southeast of the Belchers in October 2014, the average SIM inventory in the top 20 m for all stations was + 1.80 m (Table 5). This large inventory derives from sea ice that has melted through summer and, given the smaller thickness of sea ice produced locally in winter (< 1 m), likely reflects excess SIM and sea ice imported from the north by winds (Landy et al. 2017). In January 2015, the top 20 m of the water column contained an average + 0.92 m SIM (Fig. 11a; Table 5), which indicates a loss of ~ 0.88 m of SIM inventory. This change would be equivalent to the brine released by the formation of ~ 0.95 m of sea ice (correcting for the density of ice 9%). Given that the average observed ice thickness during the January sampling period was 0.5 ± 0.04 m (Fig. 11, Table 2), and deducting the ~ 12% RW in the ice implied by δ18O values (Table 7), the 0.44 m of sea ice produced locally after October can explain only about half of the brine that has been added to the top 20 m of the water column during this time period (Fig. 11a). These data imply that, like RW, brine must also have been advected into the nearshore zone southeast of the Belcher Islands between October 2014 and March 2015.

Table 7 Mean river water content (Frw) in the sea ice at various locations based on δ18O values in ice cores collected at the end of winter. The expected amount of brine (given as −SIM) that would be added to the sea surface by forming 1 m of this ice

The importance of advection is further supported by the observation that the average amount of SIM in the top 20 m further decreased by 0.5 m in the nearshore zone between January and March 2015 (Fig. 11b; Table 5). Compared to the October 2014 inventory (+ 1.8 m), the March SIM inventory reflects the loss of 1.4 m of SIM through the winter. This loss would be equivalent to brine released by the formation of 1.5 m of sea ice. Southeast of the Belchers, the average sea ice thickness in March was 0.89 m (Table 2). Correcting again for the RW content of the cores (Table 7), the average sea ice accumulation of 0.78 m between October 2014 and March 2015 is not sufficient to explain the 1.37 m of SIM removed from the water. The discrepancy is even greater if we consider the SIM inventories in the total water column (Table 6).

To the northwest of the Belchers, in January 2015, the top 20 m of the water column contained on average 1.3 m SIM (Fig. 11a, Table 5). By February–March 2015 (Fig. 11b), the SIM inventory had decreased to + 0.76 m, implying the removal of 0.54 m SIM. In January, the ice thickness was 0.50 ± 0.04 m and in March, 0.89 ± 0.05 m (Table 2), which is an average addition of 0.39 m of sea ice. Applying a fractional correction of 0.05 to account for RW in this sea ice (Table 7), yields ~ 0.37 m of sea- ice produced between January and March 2015. Although the amount of sea ice grown is about the same both northwest and southeast of the Belchers over this time period, we have no data for water column SIM in October for the northwestern waters and therefore cannot compare ice growth water inventories. Nevertheless the amount of sea ice grown in the region northwest of the Belchers after January is, once again, not sufficient to fully account for the loss of SIM inventory in the water column, which implies that the brine has been advected in from elsewhere.

Comparing SIM in waters southeast and northwest of the Belchers, we find similar inventories (not statistically different) in both January 2015 and February–March 2015, but there was a fair amount of variation between stations with, for example, SK 14 holding 1.71 m of SIM inventory (top 20 m) compared to SK 23 with 0.33 m of inventory in January (Fig. 11a, b; Table 5). This variability could be explained either by imported brine or by locally produced brine in the network of polynyas. Unlike the RW inventories, the SIM inventories do not show a statistically significant difference between the waters northwest and southeast of the Belchers. However, it is worth noting that the stations with the small SIM inventories tend to be associated with the stations with large RW inventories (e.g., SK 1, SK 21, SK 22, and SK 3).

The inventories of SIM in the top 20 m along the orthogonal transect varied widely (0.7–1.4 m) but were smaller than SIM inventories in the coastal domain in October (1.80 ± 0.2 m; Table 5), indicating brine addition between October and February–March 2015, coincident with RW addition. High fractions of SIM were present below the halocline at the stations along the transect resulting in high inventories of SIM in the total water column (3.16–4.67 m, Table 6). The SIM deeper in the water column is presumably a regional feature, remnant from summer melt.

The Relation Between River Water and Sea Ice Melt

The change in RW inventories at all stations between January 2015 and March 2015 (Table 5) is negatively correlated with the change in SIM inventories (R2 = 0.86, p < 0.001; n = 12; Fig. 13), implying that the RW imported into the region around the Belcher Islands during winter carries with it brine (negative SIM). The relationship holds true for the top 20 m inventories and for total-water column inventories but is stronger for the former. The association explains why the brine added to the water column around the Belchers during winter considerably exceeded the capacity of local ice formation to produce brine. Furthermore, the change in RW inventory and brine exhibits an approximately 1:1 relationship (slope − 1.16 m SIM/m RW with 95% CL of − 1.5 to − 0.82 and y-intercept of − 0.17 (− 0.42 to 0.07 95% CL)). We infer from the proportional relationship that the inventory of RW in the surface water sets a limit to the capacity of this water mass to hold brine. Once the brine component becomes equivalent to the RW component, further addition of brine will tend to mix and deepen the surface layer. If the brine eventually exceeds the original capacity of RW+SIM to stratify the surface layer, the water can then mix or convect more deeply. Thus, the 1:1 ratio sets a limit for the concentration of brine transported within a surface layer freshened by RW in the setting of under-ice river plumes, where wind and tidal forcings are reduced (Ingram and Larouche 1987b; Prinsenberg 1991).

Fig. 13
figure 13

Change in the sea-ice melt (SIM) inventory compared to change in the river water (RW) inventory present in the top 20 m of the study area’s water column (between January and February–March 2015)

The Source of Freshwater to the Belcher Islands Coast in Winter

The isotopic data presented here show the strong influence throughout winter of river water, especially on the southeastern side of the Belcher Islands. Furthermore, the influence of RW increases after freeze-up at some of the stations. These observations lead to two questions: what is the origin of the RW in the winter of 2014–2015, and is the increase in RW throughout winter a normal occurrence for this region? We take 2014 and 2015 data as being representative of typical winter time conditions in the southeast Hudson Bay region because freeze-up and break-up dates of sea ice in Kuujjuarapik (~ 100 km southeast of Belcher Islands) were near to the mean of the 1996–2016 period (see Fig. 9 in Andrews et al. 2018). We do not have sufficient years of data to comment meaningfully on inter-annual variability. However, it is striking that the higher RW inventories of January–February 2014 compared to January–March 2015 (see the “Distributions of RW in Winter” section) coincide with, in January–March 2014, higher sea-surface heights throughout James Bay/east Hudson Bay and stronger northerly flow (as opposed to northeasterly flow) of James Bay outflow according to the NEMO model simulation (Fig. 9). This concurrence of model results and observations provides grounds for concluding that James Bay outflow is the main control on how much river water is present in winter around the Belcher Islands and in southeast Hudson Bay. The low surface salinity (~ 25 psu) measured under the ice at Cape Jones at the mouth of James Bay throughout January–March 2015 (Fig. 6) confirms that James Bay outflow is indeed a probable source of the freshwater, and the model simulation (Fig. 9) confirms a transport pathway from the mouth of James Bay to the southern Belchers.

The salinity-δ18O data collected in our study show that surface water enriched with RW impinging on the Belchers in winter brings with it an approximately equivalent amount of brine (Fig. 13). As shown by the satellite image in Fig. 2, there was a corridor of flaw leads along the edge of landfast ice in 2015 that extended from the mouth of the La Grande River to the mouth of James Bay and then to the Belchers. This system of flaw leads would provide the setting to produce sea ice rapidly and thereby add brine to the surface layer stratified by RW input, as it transports northward along the coast and out into Hudson Bay. Thus, James Bay inflow could account for the coincident arrival in winter of RW and brine at the Belcher Islands and, therefore, account for both the increase in RW during winter and the excess of brine over the supply from local sea-ice growth. The Churchill and Nelson Rivers also supply RW in winter, but they are over 800 km away, which would require transit times in the boundary current of about 3 months, making them less likely to be the major source of the RW plume seen along the southeast coast of the Belchers in January–March, 2015.

Answering questions about longer-term change is rendered difficult due to the lack of historical oceanographic data for this region in winter. The increase of the RW component through winter runs counter to the normal seasonal hydrographs for undeveloped rivers in Hudson Bay, which see flows continuously dropping from December through to April (Déry et al. 2011). Instead, the winter increase in RW appears to match the altered seasonal hydrograph for the La Grande River as described in previous works (Déry et al. 2016; Hernandez-Henriquez et al. 2010). According to Messier (2002), monthly mean flows from the La Grande system after Phase II of the La Grande complex reach 4000 m3 s−1 (126 km3 year−1) in November and December and more than 4500 m3 s−1 (142 km3 year−1) in January and February. At flows of 3000 m3 s−1, the under-ice plume of the La Grande was found to dilute surface waters for over 100 km north of the river mouth (Ingram and Larouche 1987a). Because of reduced mixing under the landfast ice, further increases in the midwinter discharge were predicted to lead to a decrease in surface salinities and an increase in upper-layer stratification along a coastal zone some 20 km wide in northeast James Bay and extending into southeast Hudson Bay (Ingram and Larouche 1987a). The low salinity measured at Cape Jones during January–March 2015 (Fig. 6) and the low salinities and presence of river water at the southern tip of the Belchers may confirm this prediction; transit time from the mouth of the La Grande River to the Belchers (~ 220 km) estimated assuming a mean surface current of 0.05 m s−1 (Fig. 9) is 51 days, which is consistent with November–December discharge from the La Grande impinging on the Belchers in January. However, the La Grande River represents only 40% of the total river inflow to James Bay on an annual basis (~ 110 km3 year−1 out of a total 270 km3 year−1 (Déry et al. 2011)) and about 60% of the total river inflow to James Bay during the period December–March (~ 45 km3 year−1 out of a total 77 km3 year−1 (G. McCullough, personal communication)). Other large rivers such as the Moose River in western James Bay continue to flow in winter, albeit at among the lowest levels of the year (cf., Poehlman 1996). Additional tracers are needed to estimate the contribution of different rivers to the river water observed at the Belchers.

Is the strong stratification by RW observed in 2015 exceptional with respect to historical conditions? We have been able to find only one winter oceanographic dataset relevant to the Belcher Islands. In February 1977, Prinsenberg (Canadian Hydrographic Service/Canadian Centre for Inland Waters (Burlington), unpublished) collected CTD data at four stations at the southeast corner of the Belcher Islands from which it is possible to calculate total freshwater inventories and compare these with data collected in February, 2015 (Fig. 14). Excluding a very shallow station from 1977, the average freshwater inventory for the top 20 m of the water column in February 1977 was 2.8 ± 0.1 m compared to 3.4 ± 0.2 m in February 2015. Similarly, in 2015 the total water column inventories at stations of comparable depth were higher by ~ 0.5 m. While the higher freshwater inventory in 2015 is consistent with enhanced winter inflow from nearby and upstream rivers entering James Bay and Hudson Bay, the data do not provide definitive evidence of change partly because there is only the one historical year and partly because CTD data alone do not permit discrimination between RW and SIM. The SIM component for 1977 could on one hand account for the 0.5 m difference in total inventory between 1977 and 2015 or, if this water contained excess brine, could result in an even greater than 0.5 m difference in the RW inventories to balance the total.

Fig. 14
figure 14

Comparison of total freshwater inventory in the top 20 m of the water column in the coastal domain between February 1977 and February 2015. Inset shows salinity profiles

Implications for Seasonal Freshwater Cycling in Hudson Bay

The seasonal evolution of freshwater distribution in the coastal domain described above presents a different view of freshwater cycling in southeast coastal Hudson Bay than proposed by Granskog et al. 2011 (their Fig. 4.9). In their schematic diagram, the brine produced in flaw leads and polynyas located at the outer boundary of landfast ice in Hudson Bay during winter promotes convection through enhanced sea ice formation and brine production, which contributes to a deeper mixed layer and, in extreme cases, could produce water dense enough to enter the basin. An implication of our data is that the buoyancy supplied by RW in the winters of 2014 and 2015 was more than sufficient to capture the brine and transport it at the surface with the result that deep convection did not occur. Thus, the strong stratification due to RW addition may have consequences for basin deep water renewal and residence times. The emphasis in the Granskog et al. (2011) scheme was on the strength of the ice engine to produce brine, while RW and residual SIM remaining from previous summer played a smaller role in maintaining the winter surface mixed layer (WSML). Our observations for the coastal domain of southeast Hudson Bay in winter, which receives exceptional amounts of RW, suggest that the local capacity of sea ice formation to produce brine is relatively weak due to limited local sea-ice growth that is set against strong stratification supported by RW entering throughout the year.

As shown in Fig. 15, we propose that interplay between RW and brine addition upstream of the Belcher Islands determines the surface salinity and stratification observed weeks later when the RW- and brine-rich surface layer impinges on the southeast coast of the Belcher Islands. In this scheme, coastal currents transport James Bay river waters including those of the La Grande northward in a surface plume, first under the landfast ice and then northward across an extensive system of flaw leads. During transport, brine produced by rapid sea ice formation in the flaw leads modifies the fresh surface layer by increasing the density (salinity) of the water and deepening the mixed layer. However, the river water continually moving northward in the coastal current provides sufficient buoyancy to counter the negative buoyancy associated with brine addition (Fig. 15). The large winter inflow of river water in the coastal current thus prevents heat (and other constituents such as nutrients) in deeper water (> 20 m) from mixing into the surface layer, or oxygen and other properties from being delivered to deep basin waters with the sinking of dense surface waters.

Fig. 15
figure 15

Schematic of the interaction between river inflow and the formation of winter mixed layers in the flaw leads of southeastern Hudson Bay. The plan view (tilted) shows southeast Hudson Bay from the southwest of James Bay. Along the section marked A-B (red dashed line), coastal currents transport La Grande and other James Bay river water northward as a surface freshwater layer, first under the landfast ice and then northward across an extensive system of flaw leads. Brine produced by rapid sea ice formation in the flaw leads modifies the fresh surface layer as it transports northward, increasing the density (salinity) of the water and deepening the mixed layer. Within weeks, the brine-rich surface layer impinges on the southeast coast of the Belcher Islands. As illustrated in the sectional view, the large winter inflow of river water in the coastal current provides sufficient buoyancy to absorb the brine without mixing deeply, thereby preventing heat and nutrients in deeper water (> 20 m) from mixing into the surface layer, or oxygen and other properties from being delivered to deep basin waters with the sinking of dense surface waters

The increased RW supply in the southeast corner of Hudson Bay due to hydroelectric development after 1978 on James Bay rivers and possibly climate change, which, in some areas is increasing winter precipitation, would only enhance stratification. In 2015, the stratification was sufficient to trap heat below 20 m in the coastal zone (Fig. 4), but the data record for this area is far too sparse to determine whether these circumstances are new, or whether the extra RW observed in the coastal water southeast of the Belchers (compared to northwest of the islands) is contributing to changes in ice conditions, such as winter stratification of small polynyas, which permits them to freeze, or storm-driven mixing and entrainment of heat into surface waters, which could melt sea ice in winter. Further work is needed to assess whether the persistence of heat in the shallow subsurface through winter relates to strong stratification supported by RW entering throughout the year. An alternative explanation is transport from the main basin; however, there is little evidence of warming with depth in the adjacent main basin of Hudson Bay in published data. Both the winter mixed layer and bottom waters are at the freezing point (Prinsenberg 1984). The importance of even apparently minor freshening and increased stratification was also emphasized by Peralta-Ferriz and Woodgate (2015), who found that winter surface mixed layers across the Arctic have freshened and shoaled during the last 30 years, with implications for mixing nutrients and heat up into the surface layer and photic zone.

In the waters northwest of the Belchers, the seasonal evolution of freshwater distribution largely complies with the Granskog et al. (2011) conceptual model. During our study, water in this domain became well-mixed to a depth of at least 40 m (the limit of our nearshore stations) and the SIM inventory of the water column progressively decreased through the winter, consistent with brine addition. Although there was significant freshwater present at northwestern stations in winter, and episodic increases in RW at some stations (e.g., SK4), the major source of freshwater in this area was SIM remaining from the summertime, not RW. Throughout the entire study area, there is enhanced SIM inventory in late summer presumably due to the southward transport of sea-ice within Hudson Bay.

A recent review describes the dominance of stratification and small changes in the Arctic freshwater balance in determining winter mixed layer depths throughout the Arctic Ocean and the peripheral seas (Peralta-Ferriz and Woodgate 2015). A long-term assessment of mixed layer properties in Hudson Bay should be completed in view of the rather large changes to this region’s freshwater budget. The 2 years of winter data examined in this study (2014 and 2015) indicate that inter-annual differences in RW and SIM in this region also are substantial. Presumably some of the difference between years is due to conditioning in late summer and fall (i.e., the SIM and RW inventories at the start of winter) or possibly inter-annual variability in source seawater, but we do not have data for the Belchers prior to the 2013–14 winter and hence have limited our analysis here to the 2014–2015 winter. The NEMO model simulations for winters 2014 and 2015 indicate differences in both sea-surface heights and circulation of James Bay outflow around the southern Belcher Islands that require further investigation in terms of sensitivity to fall pre-conditioning, storms, or winter river discharge.

Conclusions

The origin and distribution of freshwater in southeast Hudson Bay around the Belcher Islands in winter were examined using temperature-salinity profiles and salinity-δ18O pairs, which allow differentiation of river water and sea ice melt water. Observations of coastal seawater properties in October 2014 and January–March for 2014 and 2015 were complemented by records of δ18O in sea ice cores collected in March 2015 and augmented with continuous records of the salinity just beneath the ice. Two domains exemplified by saline, unstratified waters northwest of the Belchers and fresher, stratified waters southeast of the Belchers may now be considered to be a permanent feature in Hudson Bay throughout the year. Previous studies characterized these two domains in open water, and our new data show that the two domains feature also under the sea ice during winter. The nearshore domain southeast of the Belchers is distinguished by the presence of a significant river water fraction (~ 15%). Counter to what would be expected for a system in which greatest river inflows occur during spring freshet (May–June), the inventory of river water in the water column around the Belcher Islands increases from fall through to late winter. We conclude that river water is advected into the region in winter and, using NEMO model simulations, infer that the origin of this water is James Bay. A relatively small amount of brine is produced by local sea-ice growth around the Belchers. However, brine is advected into the region along with river water and indeed a close ~ 1:1 relationship exists between accumulation of river water and accumulation of brine in surface waters during the winter. This relationship implies that the buoyancy produced by the river water is sufficient to absorb the brine produced by growing sea ice in the northern James Bay flaw lead; thus, the brine accumulates in a surface layer rather than mixing or convecting more deeply. Because of the presence of river water in the coastal domain in winter, the depth of winter convection is reduced relative to that in the offshore domain (~ 20 m vs. 40 m). With climate-change and river regulation, we anticipate increased river water flowing through and out of James Bay and southeast Hudson Bay in the winter season, which is out of cadence with the historical spring freshet. The fresh surface layer will have the capacity to store and advect brine farther away from the site of sea ice growth. The association of brine with river water is likely an important tag to define the advected river water in this region.