In this study, flux gates are set up in the FS similarly to Ricker (Ricker et al., 2018) and divided into a zonal gate (82°N, 12°W–20°E) and a meridional gate (80.5°–82°N, 20°E), as shown in Fig. 1a. To reduce the influence of SID errors, the flux gates should be set at a high latitude (Selyuzhenok et al., 2020). The monthly total sea ice volume, monthly mean air temperature, and monthly mean SST are calculated for the fixed region (Fig. 1b).

The NSIDC provides the SID data product by combining data from a wide variety of sensors, including passive microwave radiometers such as the scanning multichannel microwave radiometer, the special sensor microwave imager (SSM/I) and the SSMI sounder, and the advanced microwave scanning radiometer-Earth observing system, and International Arctic Buoy Programme buoy data. The gridded fields are on a polar stereographic projection, with a grid cell size of around 25 km × 25 km. Daily gridded ice motion data from 1979 to 2022 are obtained to calculate the sea ice volume flux, and missing SID data are replaced with data from the day before the corresponding date.

The PIOMAS is a coupled sea ice ocean model that assimilates the NSIDC near-real-time daily sea ice concentration (SIC) and daily surface atmospheric forcing and the SST of ice-free areas from the National Centers for Environmental Prediction and National Center for Atmospheric Research reanalysis. The original PIOMAS data are gridded on a 25-km EASE-2 grid to maintain the spatial resolution of the SID data.

Monthly upper 300 m ocean heat content (OHC) is available at the ocean reanalysis system 5 (ORAS5) provided by the European Centre for Medium-Range Weather Forecasts, covering 1979 until the present. Monthly 2 m SAT data and monthly SST data are obtained from ERA5 to analyze their effects on the sea ice in the GS.

The large-scale atmospheric circulation indexes used in this study include AO [the first mode of the empirical orthogonal function (EOF) of north of 20°N sea level pressure (SLP)], NAO (the first mode of the EOF of North Atlantic SLP), and DA (the second mode of EOF of SLP within the Arctic Circle north of 70°N). Monthly AO and NAO can be obtained from the National Oceanic and Atmospheric Administration, while DA is provided by Bi Haibo (personal communication).

The volume flux through FS is calculated, and the gates are divided into a zonal gate (82°N, 12°W–20°E) and a meridional gate (80.5°–82°N, 20°E). The fluxs at both gates (FV_lon and FV_lat) are calculated as Eqs (1) and (2), then summed to get the total volume flux (FV_fram) through FS as Eq. (3).

where $ {V}_{i} $ is the SID speed, $ {T}_{i} $ is the effective SIT, $ \Delta x $ is the grid width, and $ {\mathrm{\mathit{N}}_\mathrm{m}} $ and ${\mathrm{\mathit{N}}_\mathrm{z}} $ refer to the number of grids in the meridional and zonal gates.

The uncertainties of SID and SIT products determine the error in sea ice volume flux. The uncertainty of SID is stored in each grid of data, while that of SIT is estimated to be about 0.19 m (Zhang and Rothrock, 2003). The uncertainties of the daily field and annual volume flux (AVF) are estimated as follows (Kwok, 2009):

where $ {{\sigma }}_{\mathrm{d}} $ is the uncertainty of the daily SIM data, $ {{\sigma }}_{\mathrm{F}} $ is the daily field uncertainty, $ L $ is the length of the flux gate, N is the number of independent samples across the gate, and $ {N}_{D} $ is the days for the period examined. The average volume flux uncertainties through FS are estimated to be 6.42 km³ during winter and 3.87 km³ during summer, which are less than 10% of the seasonal magnitude, and thus, can be negligible.

The sea ice volume in the GS (GV) is calculated by multiplying the effective SIT by the grid area, as shown in Eq. (6).

where N_GS refers to the number of grids in the GS.

The monthly sea ice volume variations in the GS are obtained by subtracting the sea ice volume in the current month from that in the following month [Eq. (7)]. Advection process and local processes can determine sea ice volume variations [Eq. (8)]. The sea ice volume variations for advection can be expressed by the monthly volume flux (MVF) through FS as Eq. (9), and the local processes can be obtained by subtracting the advection from the total sea ice volume variations each month, which can be calculated as Eq. (10).

where $ \Delta V $ denotes the monthly change in the sea ice volume in the GS, $ {V}_{{\mathrm{a}}} $ denotes the change in the sea ice volume owing to sea ice inflow through the FS in the north (a positive value represents sea ice from the north to the south, and a negative value represents its northward movement), $ {V}_{{\mathrm{l}}} $ denote the change in the sea ice volume owing to the local processes (a positive value represents an increase in the sea ice in the GS, and a negative value represents a decrease in the sea ice in the region), $ {V}_{m} $ denotes the sea ice volume in the GS for the month, and $ {V}_{m+1} $ denotes that in the subsequent month, and D denotes the number of days in the month. The sea ice volume flux through the Denmark Strait is less than 2% of that through FS, and thus, the flux of the southern boundary is not considered in this study.

Compared with the sea ice volume flux obtained by other researches, our findings exhibit a favorable agreement with the observed trend (Fig. 2), and the obtained correlation coefficients surpass 0.69 (Table 1). However, our values appear to be generally lower than those of other studies. This discrepancy can be attributed to variations in the data sources. Specifically, the PIOMAS model exhibits systematic differences by overestimating the thickness of thin ice and underestimating the thickness of thick ice and the GS dominated by thick ice (Schweiger et al., 2011). In addition, the SID velocity derived from passive microwave data may be slower than that obtained from SAR data, which may result in low flux results. Moreover, the positioning of gates in the FS varies in different studies. Generally, the farther south the gate position, the larger the error in the thickness and in the drift rate (Kwok et al., 2004). The use of different algorithms can also introduce biases in the calculation results.

Owing to differences in the data sources and in the research periods, the volume flux through FS calculated by different scholars demonstrates different trends (Kwok et al., 2004; Selyuzhenok et al., 2020; Li, 2021; Wang et al., 2022). Some scholars calculated the volume fluxes through FS during the period of 1978–2018 and concluded that they demonstrated no significant trend (Li, 2021). Meanwhile, some scholars calculated the volume fluxes from 1979 to 2017 and observed an increasing trend (Selyuzhenok et al., 2020), whereas others calculated the volume fluxes during the period of 2010–2019 and discovered a decreasing trend (Wang et al., 2022). Owing to the lack of reliable SIT data, the volume fluxes obtained by different scholars exhibit large differences, making obtaining a consistent trend difficult (Bi et al., 2018). On the one hand, Arctic sea ice is melting at an accelerated rate (Spreen et al., 2020). On the other hand, the continuous decrease in the SIT (Rampal et al., 2009; Spreen et al., 2011; Howell et al., 2013) has contributed to the SID accelerating trend. Therefore, the decrease in the SIT and the increase in the SID can simultaneously affect the sea ice volume flux; however, which of the two parameters is more decisive has yet to be determined.

From 1979 to 2022, the annual mean sea ice volume flux through FS was 136 km³/month. As depicted in Fig. 3, the highest flux occurred in 1995, with an annual mean volume flux of 231 km³/month, and the lowest flux occurred in 1984, with an annual mean volume flux of only 48 km³/month. The volume flux through FS exhibited considerable variability before 2000. Notably, anomalously high volume fluxes occurred in 1989 and 1995, which aligns with the observations of previous studies (Arfeuille et al., 2000; Lindsay and Zhang, 2005). This study suggests that the anomaly in 1989 was due to an increase in the SIT, whereas the anomaly in 1995 was caused by the accelerated southward drift of the sea ice (Lindsay and Zhang, 2005). The volume flux stabilized after 2000, but two instances of very low fluxes were recorded in 2016 and in 2018. However, no significant trend is observed during the two periods. Sumata et al. (2022) found that the minimum output in 2018 should be attributed to the decline of SIT as well as regional scale atmospheric anomalies. The drop of effective SIT was 0.74 m, corresponding to a reduction of sea ice volume by 55% in FS. Additionally, the reversals of zonal SLP gradient across FS led to southerly winds that further reduced the annual mean southward motion of sea ice to FS. Similar to 2018, the annual mean effective thickness in 2016 was only 0.72 m, which is significantly thinner than the average in 2003–2018, resulting the lower export (Sumata et al., 2022).

The sea ice volume flux through FS shows significant seasonal variations, with higher fluxes observed during winter (Fig. 4). In this study, winter is defined as the period from October to April, and summer refers to the months from May to September. On average, the flux in summer contributing 27.37% annual total while the winter accounting for 72.63%.

As shown in Fig. 5, the sea ice volume in the GS from 1979 to 2022 exhibites no significant trend. However, a notable downward trend emerged during the period of 2001–2022, which is characterized by a decline rate of −5.36 km³/a. Throughout the entire period of 1979–2022, all the months, except June (downward trend but do not pass the significance test), show a significant downward trend. Notably, the winter months experienced a highly pronounced and rapid decrease, with January and March exhibiting a decreasing rate of −4.19 km³/a and −4.46 km³/a, respectively. Whether on an annual or on a monthly scale, the sea ice volume in the GS declined rapidly during the second period. Some studies demonstrated that the SIT has more influence on the sea ice volume compared with other factors, particularly during winter, when it can account for over 80% of the variability (Li, 2021). Thus, the accelerated decrease during the period of 2001–2022 indicates that the sea ice in the GS has been undergoing rapid thinning since the beginning of the 21st century (Liu et al., 2020; Li, 2021).

The changes in the sea ice volume in the GS are influenced by a combination of advection and local factors. The sea ice volume in the GS exhibits a decreasing trend, whereas no significant trend can be seen in the sea ice volume flux, which indicates that sea ice melting is increasing rapidly. When the sea ice flux replenishment fails to offset the sea ice melting, it may reduce the sea ice volume, which shows that the local processes on sea ice melting in the GS is strengthening.

As depicted in Fig. 6a, the distribution of the effective SIT in the GS exhibits a V-shaped pattern from north to south, with thick sea ice concentrated primarily in the northern region. During the period of 1979–2000, the sea ice was distributed extensively, and the effective SIT ranged of (1.00 ± 0.50) m (Fig. 6b). However, compared with the first period, the sea ice extent decreased noticeably and shifted northward during the period of 2001–2022 (Fig. 6c). In addition, the effective SIT decreased significantly within the range of (0.05 ± 0.03) m.

During the period of 1979–2022, the rate of change in the effective SIT in the GS showed a decreasing trend, with a T-distribution (Fig. 7a). During the period of 1979–2000, the areas that exhibited thickness changed rapidly, primarily around 12°W, within the range of (−0.01 ± 0.02) m/a (Fig. 7b). The spatial distribution of the decreasing SIT trend during the period of 2000–2022 was also a T-distribution, ranging of (−0.02 ± 0.04) m/a (Fig. 7c).

During the period of 1979–2000, the SST and SAT exhibited an upward trend but do not pass the significance test. However, during the period of 2001–2022, the SST and the SAT showed a significant increase. The years with extreme values generally coincided with one another, with maximum values in 1984, 1990, and 2016 and minimum values in 1988 and 1997. The warming of the SST and the SAT at the beginning of the 21st century accelerated sea ice melting, which corresponds to the accelerated decline of sea ice observed in the second period.

As depicted in Figs 8b and c, the SAT reached the highest value in July and the lowest value in February, whereas the SST generally reached the highest value in August and the lowest value in March. During the months with extreme temperatures, a one-month lag can be seen in the SST compared with the SAT. This lag can be attributed to the larger specific heat capacity of the ocean compared with that of the atmosphere (Li, 2021). The ocean absorbs the same amount of energy but undergoes a relatively small temperature change. In addition, temperature changes from water mixing in the ocean require a long period of time. The slow heat transfer process caused by ocean currents can contribute to the delayed response of the SST to external changes. The combined factors can result in a one-month lag in the SST compared with the SAT.

As depicted in Fig. 9a, the spatial changes in the SST from 1979 to 2000 show a decreasing trend in the northwest region and an increasing trend in the southeast region, with the most noticeable warming areas located along the coast of Svalbard and in the Denmark Strait. During the period of 1979–2000, the changes in the SST were relatively small, with a warming area located on the southern side of Svalbard and few cooling areas (Fig. 9b). However, the changes in the SST were highly pronounced during the period of 2001–2022 (Fig. 9c). The spatial distribution exhibits cooling in the northwest and warming in the southeast. Notably, the warming area shifted significantly northward in the second period compared with that in the first period, and the rate of warming demonstrated a significant increase.

The spatial distribution of the SAT changes appears to be highly consistent and smooth (Figs 9a–c). During the different periods, a general warming trend can be observed in the northeast and southwest regions, whereas a cooling trend is evident in the northwest and southeast regions.

To facilitate the comprehensive comparison of the impact of ice advection from the north and local processes on sea ice changes, this study applies a sign convention to the volume variations attributed to local forcing. Specifically, in Fig. 10, the positive values indicate a reduction in the sea ice volume in the GS, whereas the negative values represent an increase in the volume.

Sea ice volume variations influenced by advection and local forcing exhibit similar trends (Figs 10a–c). During winter, the impact of ice advection from the north on sea ice volume variations is typically greater than that of local processes (Fig. 10b). Conversely, during summer, the local processes have a greater influence (Fig. 10c). However, throughout the different periods, the sea ice advection primarily replenishes the sea ice, whereas the local forcing predominantly melts the sea ice. The balance between sea ice volume export to the GS through FS and regional changes is strongly influenced by the sea ice transport and local processes. The volume of sea ice in the GS and sea ice flux are not strongly correlated, indicating that the sea ice variability in the GS can be greatly affected by the local processes (Selyuzhenok et al., 2020; Chatterjee et al., 2021).

The variations in the sea ice volume induced by local factors are strongly linked with the local sea ice dynamics process and oceanic conditions (Campbell et al., 1987; Visbeck et al., 1995; Kern et al., 2010; Selyuzhenok et al., 2020; Chatterjee et al., 2021). The sea ice volume variations caused by ice advection from the north do not exhibit a significant trend, whereas the total sea ice volume in the GS displays a significant decreasing trend. Therefore, the influence of local factors may be masked by the overall decreasing trend of the sea ice volume. A highly accurate assessment of the impact of local factors on sea ice can be conducted by considering the rate of sea ice changes induced by local factors to the total sea ice volume. As shown in Fig. 11, before 2000, the influence of local factors on sea ice was relatively limited, with the rate remaining below 1. However, since the start of the 21st century, the influence rate has exhibited a sharp increase, reaching 2.30 in 2012 and 2.73 in 2016, which indicates the intensified local processes. The average annual Arctic sea ice volume in the GS in 2016 and 2012 is the minimum and sub-minimum (Fig. 5), respectively. The loss of sea ice led to more open water, along with anomaly high SAT and SST in 2012 and 2016 (Fig. 8), accelerating the melting of sea ice. As sea ice melts, the open water area in the GS will expand, which will lead to increased energy absorption by the ocean and further melting of the sea ice. On the other hand, low SLP over GS associated with stronger cyclonic Greenland Sea Gyre (GSG) and stronger Ekman drift, which will weaken the upper-ocean stratification and enhanced vertical mixing of anomalous warm Atlantic Water (AW) by reducing the freshwater in the GS. Under the weak stratification condition, the AW will warm the surface waters by enhanced vertical mixing and then further inhibiting new sea ice formation or melting the sea ice from bottom (Ivanova et al., 2012). The increased salinity caused by the vertical mixing further promotes a stronger GSG circulation and then reduce sea ice in the region, completing a positive feedback mechanism (Lauvset et al., 2018; Chatterjee et al., 2021).

The sea ice volume in the GS exhibits a strong correlation with the SAT and the SST. The SAT and the SST demonstrate a significant increasing trend (Fig. 12) and have a negative correlation with the sea ice volume in the GS.

As depicted in Fig. 12, the impact of the SAT on the sea ice volume exhibits a lag of approximately 2 months to 3 months, and that of the SST exhibits a lag of 1 month to 2 months. In terms of lag, the SAT exhibits the earliest change, thereby influencing the SST, which in turn indirectly affects the sea ice through the SST. The sea ice growth and melting can be determined by the heat exchange processes occurring within, on top, and at the bottom of the sea ice. Such processes are jointly controlled by the SAT and the SST. When the SST reaches the freezing point (−1.8℃), the density of the sea will increase, then decrease, and warm water from the Atlantic Ocean will rise to the surface and exchange heat with the surface air once again. Sea ice will not begin to melt until the warm water that rose to the surface reaches the freezing point (Li, 2021). In addition, sea ice exchanges heat with the environment as it grows and melts, which in turn can affect the changes in the SAT and in the SST.

The SST and SAT anomalies are analyzed using the EOF to obtain the spatial distribution and corresponding time series coefficients of their first two modes. The first mode of the SST captures the interdecadal scale changes, and when combined with the time series coefficients (Fig. 13b), it reveals two contrasting anomalies. Specifically, anomalously low SST can be observed during the period of 1979–2000, followed by anomalously high SST during the period of 2000–2022. The areas experiencing rapid warming are concentrated along the coast of Svalbard and on the eastern side of the Denmark Strait (Fig. 13a), which correspond to the areas of the Irminger current and the West Spitsbergen current. Conversely, the slow warming observed along the eastern coast of GS can be attributed to the presence of the East Greenland current, which transports cold water from the Arctic. The first mode of the SAT exhibits characteristics similar to those of the SST (Fig. 13e), with anomalously low temperatures before 2000 and anomalously high temperatures after 2000. However, its interdecadal characteristics are less pronounced than those of the SST (Fig. 13f), because the SAT is more susceptible to climatic events, which may result in many frequent fluctuations.

The second mode of the SST captures the interannual variability in the SST and reveals that the regions with high SST values and positive time coefficients align closely with the areas with a rapidly decreasing effective SIT (Fig. 7a). The regions form distinct stripes that extend southward from the northern part of the Svalbard Archipelago (Fig. 13c). The least squares fit analysis of the time coefficients of the second mode after 2000 (Fig. 13d) reveals a significant upward trend (0.04). This finding indicates that the warming pattern observed in the striped region occurred frequently after 2000, which corresponds to the rapid decline in the sea ice volume in the GS during the period. When combined with the time coefficients (Fig. 13h), it shows that the northwest region experienced anomalous warming during the period of 1979–2000, whereas the southeastern region experienced anomalous warming during the period of 2001–2022. The second mode of the SST exhibits a spatial distribution characterized by a seesaw pattern of opposite changes in the northwest and southeast regions (Fig. 11g), which has been more variable after 2000. The oceanic temperature may be responsible for the sea ice variation in the GS, which consistent with other study that show the sea ice volume loss in the GS are largely determined by the oceanic conditions (Selyuzhenok et al., 2020).

The AO index has a considerable influence on the volume flux. Numerous studies highlighted the impact of the AO on sea ice motion and emphasized its significant association with sea ice output from FS (Zhang et al., 2003, 2023; Kwok et al., 2013; Howell and Brady, 2019). During a positive AO, the sea surface pressure in the Arctic region decreases, which will weaken the Beaufort anticyclonic circulation. This weakening will enhance the cyclonic tendency of the sea ice movement, resulting in a westward shift of the transpolar drift centerline. Consequently, an increased amount of sea ice, particularly from the northern region of the Canadian Arctic Archipelago, will flow out through FS (Rigor et al., 2002). Conversely, during a negative AO, the intensification of the Beaufort anticyclone will facilitate the transport of sea ice from the western part of the Arctic Basin to the east through enhanced easterly anomalous winds. This process can suppress the export of sea ice (Rigor and Wallace, 2004). A significant correlation coefficient (r = 0.47) can be observed between AO and the FS volume flux before 2000 (Fig. 14a). However, this relationship became insignificant (r = 0.10) after 2000 (Fig. 14b). The correlation analysis reveals a weakening relationship between the AO and sea ice output during winter. The reduction in sea ice thickness may be responsible for the decreased correlation, as it alters the dynamics of sea ice and making it less responsive to AO variations. The effect of NAO on sea ice output is similar to that of AO but less significant.

DA is also an important driving force for the transport of Arctic sea ice to the North Atlantic. When DA is in positive phase (positive SLP anomaly in the Canadian Islands, negative anomaly in the Barents Sea), the meridional wind blows from the west to the east and enhances the transpolar drift, increasing Arctic sea ice output (Wu et al., 2006; Wang et al., 2009; Zhang et al., 2022). In terms of correlation (Fig. 14c), DA exerts greater effect on sea ice output than AO (NAO). The positive phases of DA consistently promote sea ice export during summer (r = 0.37) and winter (r = 0.43) in 1979–2022.

Since SST is closely correlated with the local sea ice in the GS (Section 4.2), we calculated the correlation coefficient between SST and AO/NAO/DA index. The results indicate that NAO has greater impact on SST (r = 0.26), while AO (r = 0.01) and DA (r = −0.13) are not significant. NAO is the leading mode of the pressure variability in the northern North Atlantic (Woollings et al., 2010; Moore and Renfrew, 2012; Foukal and Lozier, 2017; Selyuzhenok et al., 2020), and plays a crucial role in regulating ocean circulation in the NS. Specifically, during the positive phase of the NAO, the cyclonic atmospheric circulation over NS become stronger and results in the cyclonic oceanic circulation in NS intensifies (Schlichtholz and Houssais, 2011). This leads to weaker upper-ocean stratification and stronger vertical mixing of anomalous warm AW, increasing the temperature in the upper ocean (Section 4.1).

To investigate the influence of NAO on the ocean circulation, we calculated upper 300 m OHC over the GS. The reduction of sea ice in the GS are partially driven by the increasing of OHC. As illustrated in Fig. 15, OHC exhibits no significant trend in 1979–2000. After 2000, however, an increasing trend in OHC was highly pronounced coincides with the loss of sea ice in the GS. Compared to the period before 2000, the correlation coefficient between winter NAO and OHC increased from 0.08 to 0.40 during the period of 2001–2022 (Table 2), which indicated that the winter NAO has a greater influence on the oceanic systems after 2000. This finding is consistent with previous studies revealed that a stronger oceanic heat advection during the positive NAO (Walczowski, 2010; Raj et al., 2018). This suggests that the positive winter NAO may accelerate the melting of sea ice by increasing OHC. This is confirmed by Ivanova et al. (2012), who exhibited enhanced bottom melting in this region associated with positive NAO.