Increasing fluctuations in the Arctic summer sea ice cover are expected with future global warming

The loss of Arctic sea ice (ASI) represents a major transformation in the Arctic region, impacting regional and global climate, ecosystems, and socio-economic structures. Observational and reanalysis data have consistently shown a notable shift in polar environmental conditions over recent decades, marked by a substantial reduction in the ASI area and a rise in the variability in its coverage and distribution. Utilizing data from the latest Coupled Model Intercomparison Project phase, our study reveals a consistent pattern highlighting a fundamental shift in ASI dynamics preceding total loss. We observe increasing fluctuations in the September ASI area as the threshold for an ice-free Arctic is approached across various scenarios and models. This pattern is particularly concentrated in the Central Arctic (CA) sub-region. Spatial analyses reveal increasing variance along the CA’s northern coastlines, accompanied by a substantial increase in open water coverage, underscoring the shift from stable to highly variable ice conditions in this region. Additionally, our findings suggest a potential link between increased ASI fluctuations and variability in surface wind speeds. These specific results underscore the urgency of multidisciplinary approaches in addressing the challenges posed by ASI variability, with implications for marine ecosystems, Indigenous communities, and navigational safety.


Introduction
The Arctic is a critical component of the Earth's climatic system and experiences rapid environmental changes with far-reaching implications.Characterized by its extensive sea ice cover, this region plays a pivotal role in influencing climatic phenomena occurring both locally and remotely [1][2][3][4].
The Arctic sea ice (ASI) has undergone substantial transformations in recent decades, a trend that represents one of the most notable environmental shifts in the Arctic.Since the late 1900s, there has been a consistent and well-documented decline in ASI across every month and region, with considerable changes in its composition and thickness [5][6][7][8][9][10].The proportion of multi-year, thicker ice has dramatically decreased, while first-year, thinner ice now dominates the ASI landscape, making the Arctic increasingly susceptible to rapid changes and more frequent seasonally ice-free conditions [11][12][13].Recent observational and reanalysis data indicate shifts in spatial patterns and seasonal behaviors, highlighting an ongoing and considerable decrease in ASI along with a rise in the variability of its extent and distribution [14][15][16].
The loss of ASI is a key driver of Arctic amplification [17].This phenomenon, where the Arctic warms more rapidly than the global average, is well-documented in both climate models and observational data [18][19][20].Arctic amplification is driven by various feedback mechanisms, with the sea ice-albedo and lapse-rate feedbacks playing particularly important roles [21][22][23][24].The sea ice-albedo feedback intensifies warming as rising surface temperatures lead to the melting of reflective ice and snow, causing the darker ocean water to absorb more heat and further accelerate ice melt.The lapse-rate feedback involves uneven atmospheric heating, with the Arctic experiencing less heat loss due to stable air layers near the surface, unlike the tropics where heat is more efficiently emitted to higher altitudes.Arctic warming is further accelerated by increased long-wave radiation and heat fluxes over open waters during colder periods [17].
Currently, the Arctic is undergoing a transition towards ice-free summers, with the possibility of such occurrences happening before the middle of this century [5,25].A shift towards an ice-free Arctic summer would have profound impacts on both regional and global climate and ecosystems, threatening the stability of northern Indigenous communities and wildlife [3,[26][27][28][29][30][31][32][33].Additionally, the disappearance of the ASI opens up new possibilities for shipping routes and resource exploitation.While this shift might be considered beneficial from an economic perspective, it could lead to significant social and environmental challenges [34][35][36][37].
Given the impacts of an ice-free summer Arctic, it is crucial to project and understand possible future trajectories of this critical region under different climatic scenarios.Here, we analyze September ASI area data from the latest phase of the Coupled Model Intercomparison Project (CMIP6) to investigate future changes in the summer ASI.We focus on four different Shared Socioeconomic Pathways (SSPs) future scenarios and combine them with the corresponding historical runs [38].The four considered future scenarios correspond to very high (SSP5-8.5),high (SSP3-7.0),intermediate (SSP2-4.5),and low (SSP1-2.6)levels of greenhouse gas emissions.Each of these pathways presents a different probability of exceeding the 1.5 • C and 2 • C global warming levels relative to pre-industrial temperatures (1850-1900) [39].In the following, we will refer to the combined historical and future scenarios by using only the names of the SSPs.
Despite the notable inter-model spread [25], our analysis reveals a consistent pattern among the majority of models within the CMIP6 ensemble, both on a Pan-Arctic and a local scale.Specifically, these models demonstrate an increase in the variability in the ASI cover during summer as it converges towards an ice-free state, defined as an ASI area below 10 6 km 2 .

Methods
In our study, we employ ASI area data from the CMIP6 ensemble.For the total ASI area, we utilize either the variable for Northern Hemisphere sea ice area siarean or we multiply the variable for sea ice concentration siconc times the individual grid cell area areacello and then sum over the Northern Hemisphere.For the spatial analysis, we use the siconc variable.We analyze models that have at least one available run for each of the future scenarios considered.We do not exceed three runs per model to avoid biases that might arise from over-representing data from any single model when averaging over the ensemble.We combine the SSP runs with the corresponding historical ones, to obtain 251-year-long time series from 1850 to 2100.We then exclude the experiments for which the threshold for an ice-free Arctic (10 6 km 2 ) is not reached, as well as the runs that reach the ice-free state in the past (table S1).

Year of complete sea ice loss
In some model experiments, the ASI area fluctuates around the 10 6 km 2 threshold for a few years after the sea ice has dropped below this level for the first time in September.Therefore, we smooth the time series using a 20-year centered moving average and determine the year of the complete ASI loss from the smoothed curve (table S2).

Time synchronization
To consistently derive multi-model ensemble averages, we implement a standardized method each time we seek to obtain an ensemble mean.This involves gathering all the collected model experiments and aligning them temporally based on their projected year of complete sea ice loss (as shown in figures S1-S4 for the ASI area time series).By adjusting the timelines of the models to synchronize in this critical year, we ensure a uniform reference point for all simulations.Subsequently, we calculate the average of the time-aligned curves to produce a comprehensive, multi-model ensemble representation (figures 1(a) and S5-S7).We apply the same method when computing multi-run means for single models.

Variance
We apply a 20-year centered moving-average smoothing to obtain the residuals for the ASI area for each experiment, that is, the difference between the ASI area data and the moving-average smoothing.The variance of the residuals is estimated using windows of different lengths.We apply the same time synchronization method used for the ASI area curves and then average the variance curves to obtain the multi-model ensemble mean.Results for the variance of the residuals within a 50-year-long window under ).In the supplementary material, we show results for 75 and 100-year-long variance windows and for all the remaining scenarios (figures S8-S11).The multi-run variance-mean for each model is computed analogously.
To understand the dynamics leading up to the complete loss of the ASI, our subsequent analysis will concentrate specifically on the 50 years leading up to the ice-free state.
To accurately evaluate the trends within the datasets, we employ a least squares approach to fit a linear trend to both the ensemble mean and the multi-run mean-variance curves, subsequently calculating the slope of each fit (tables S3-S7).The regional analysis of the multi-sensor analyzed sea ice extent (MASIE) sub-regions of the Arctic is performed analogously [40] (figure 2).
The statistical significance of the trends is determined by computing 10 000 Fourier surrogates with random phases of the variance for each scenario, both averaged over the ensemble and for each model.We then derive p-values for the slopes measured in the model data by computing the percentage of surrogates with a slope larger or equal to the original data slope and consider the increase significant when p < 0.05.We use the same significance test for the MASIE regions analysis.
For the grid cell analysis, we compute the variance of the ASI area in each grid cell for each model run and each scenario.We show results for the variance computed on a 50-year-long window (figures 3 and S12-S14).To emphasize the sub-regions with the greatest variability, we normalize the variance across the entire Arctic region at each time step.This approach allows us to focus on identifying areas of high variability without concerning ourselves with the specific values of the variance.We subsequently synchronize the runs from each model and average them.To obtain an ensemble mean, we regrid each model on a 325 × 325 Lambert-Azimuthal equal-area grid using the climate data operators, then synchronize and average [41].At each time step, we mask out the grid cells whose sea ice concentration drops below 15%.

Open water
The open water area in the Central Arctic (CA) sub-region is computed as the area of the grid cells whose sea ice concentration falls below 15%.This computation is performed for each run under each scenario.We subsequently synchronize the results from the runs from each model and average to obtain a multi-run curve for each model.Finally, we express the open water areas as a percentage, representing the ratio of the open water surface area to the total area of the CA (tables S8-S11).

Probabilities
To estimate the sea ice dynamics in the CA, we calculate the probability of sea ice absence (sea ice concentration below 15%) within individual grid cells for each model.This calculation is conducted over 5-year intervals throughout the fifty years leading to the complete ASI loss.We then quantify the extent of  these ice-free regions, expressing them as a percentage of the total CA's surface area.To obtain an ensemble mean, we remap each model field to a 325 × 325 Lambert-Azimuthal equal-area grid, then synchronize and average [41] (figures 4 and S15-S17).Percentage of the total CA's surface covered by grid cells with different probabilities of being ice-free (defined as sea ice concentration below 15%) over 5-year intervals, averaged over the model ensemble.Fifty years before the total ASI loss, the majority of the CA's surface comprises grid cells having a 0 probability of becoming ice-free in the subsequent five years.After a few decades, we observe a transition to a new, less stable state where most areas are characterized by grid cells having a probability equal to or larger than 0.5 of being ice-free.

Surface wind speed
We compute the wind speed for the model runs that reach the threshold for an ice-free Arctic in each scenario using the near-surface wind variables, namely uas for eastward winds and vas for northward winds.We consider an average of the summer months of July, August, and September, and we restrict to latitudes above 60 degrees north.We detrend and calculate the variance of the residuals over each grid cell using the same window lengths employed for the ASI area for each separate run.We subsequently synchronize and average to obtain multi-run means for each model.For each scenario, we categorize the models into two distinct groups based on whether they exhibit increased fluctuations in the ASI area or not (tables S4-S7).After remapping to a 325 × 325 Lambert-Azimuthal equal-area grid, we compute the average surface wind speed variance among the models in each group at various time steps.We then subtract the results of one group from the other to highlight differences (figure 5).The statistical significance of the observed trends is evaluated through a permutation-based significance test.This test is conducted to ascertain if the noted differences in variance between the two groups result from our specific choice of model grouping or if similar differences would emerge under alternate configurations, as well.We generate 10 000 different model grouping configurations for each time step and window length.We then compute the sum of the squared differences (SSD) between the two groups by summing over each grid cell i as: The p-value for the SSD measured with the original model grouping is derived by computing the percentage of grouping configurations with an SSD larger or equal to the original one.We consider the difference significant when p < 0.05 (figure S18).

Results
As the complete loss of the ASI is approached, we observe increasingly large fluctuations in the ASI area in the CMIP6 ensemble experiments (figures 1(b), S8-S11 and table S3).The majority of the models show increasing variance in the ASI area as they reach the ice-free state in each future scenario (69% of the models in the SSP5-8.5, 65% in the SSP3-7.0,64% in the SSP2-4.5, and 79% in the SSP1-2.6).This feature is robust with respect to varying choices of the window length chosen to compute the variance (tables S4-S7).
To better understand the spatial distribution of this phenomenon, we repeat the analysis on 10 different MASIE sub-regions.The regions we consider are the Beaufort Sea, the Chukchi Sea, the East Siberian Sea, the Laptev Sea, the Kara Sea, the Barents Sea, the Greenland Sea, the Baffin Bay, the Canadian Archipelago, and finally the CA (figure 2).We exclude the model AWI-CM-1-1-MR from this regional analysis due to the unstructured nature of its grid and the errors in the ASI area resulting from the remapping process.Across all scenarios and for all window lengths, the CA consistently shows increasing variance as the ice-free state is approached.In particular, it occurs in 96% of the models in the SSP5-8.5 and SSP3-7.0future scenarios, and 100% of the models in the SSP2-4.5 and SSP1-2.6 scenarios.While we also find increasing fluctuations in other regions, this phenomenon is observed in less than half of the models examined (table S12).
To quantify the magnitude of the variability in the ASI area on a refined local scale, we measure the variance for each grid cell.We find a notable shift in the region exhibiting the highest variance as the total loss of the ASI is approaching.Initially, this region encompasses the coastal areas of the North American and European continents.However, as the ASI loss progresses, the region of highest variance transitions to the northern coast of Greenland and the region above the Canadian Archipelago (figure 3).Both of these areas are part of the CA sub-region in the MASIE classification (figure 2).
As the fluctuations in the CA sub-region increase, we observe a profound and rapid transformation in this region.Over the five decades leading to the complete ASI loss, the models predict a pronounced reduction in ice coverage, transitioning to a new state where open water now constitutes the majority of the CA's surface.To quantify this transition, we analyze the mean percentage of open water relative to the total area of the CA.Averaged over the ensemble, this reveals a clear increase in open water coverage over the 50 years preceding the total ASI loss.The mean open water coverage in the period that spans from 100 to 50 years before the complete ASI loss is 2.4% under the SSP5-8.5 scenario, 2.9% under the SSP3-7.0scenario, 3.8% under SSP2-4.5, and 4% under SSP1-2.6.However, in the year of complete loss of the ASI we observe a substantial shift.The fraction of open water increases to 63%, 68%, 65%, and 64% of the CA's area, under each respective scenario (tables S8-S11).These percentages are in stark contrast to the historical averages, underscoring the profound and rapid environmental changes underway in the Arctic.
The combination of increased open water coverage and fluctuations in the sea ice area of the CA sub-region contributes to reduced predictability of its state.For each model, we calculate the probability of sea ice absence (defined as sea ice concentration below 15%) in individual grid cells over 5-year intervals during the five decades preceding the complete ASI loss.We then quantify the area of these regions with absent sea ice as a percentage of the total CA's surface.Averaging across the ensemble, a distinct shift in the distribution of these areas is evident.Fifty years before the total ASI loss, a relatively stable state is observed.In this period, the majority of the grid cells in the CA have a 0 probability of becoming ice-free in the subsequent five years (82% under the SSP5-8.5 scenario, 77% under the SSP3-7.0scenario, 71.1% under SSP2-4.5, and 68.3% under SSP1-2.6).Within a few decades, this state transitions to one where most areas are characterized by grid cells having a probability equal to or larger than 0.5 of being ice-free (figures 4 and S15-S17).In particular, in the year of complete ASI loss, most of the CA's surface exhibits very high variability.Specifically, 61.1%, 62.8%, 62%, and 70.1% of the CA, respective to each scenario, is covered by grid cells with a probability equal to or larger than 0.5 of being ice-free.
We explore the potential association between the increasing fluctuations in the ASI area and the variations in surface wind speeds within the latitudes encompassing the polar cell (above 60 degrees north).We consider the mean surface wind speed for the summer months of July, August, and September, and assess the variance across each grid cell using the same window lengths previously employed for the ASI area analysis.Notably, the models displaying increased variance in the September ASI area manifest a generally higher mean variance in surface wind speed compared to those models that do not exhibit increased variability in the ASI under all evaluated scenarios at different time steps, especially over the CA region (figure 5).Using permutation-based significance tests, we find that this pattern is statistically significant (p-value < 0.05), for the SSP5-8.5 and SSP3-7.0scenarios (figure S18).However, the significance of this trend does not extend to the other scenarios (see Methods section for details).We also conduct the analysis with a focus solely on the ocean surface, deliberately excluding continents from consideration.The results and statistical significance are analogous.

Discussion
We analyze CMIP6 projections for four different SSP scenarios, combined with the corresponding historical runs.Despite the large inter-model spread, varying sea ice modeling approaches, and the different forcing levels represented by the pathways, our results reveal a consistent pattern across the model ensemble.As the complete loss of the ASI approaches, we find increasing fluctuations in the ASI area, especially in the CA sub-region (figures 1(b), S8-S11, tables S3 and S12).These findings are robust across various models and scenarios, suggesting a wide and fundamental shift in ASI dynamics.
While we cannot conclusively establish a direct causative link between changes in surface wind speeds and variations in the ASI area, our analysis reveals a possible link between the increasing variance of surface wind speeds and the increasing variance of the ASI area (figure 5).The observed relationship between increased variance in surface wind speeds and ASI area fluctuations, particularly under the SSP5-8.5 and SSP3-7.0scenarios (figure S18), suggests that atmospheric dynamics play a significant role in the Arctic's climate system.The lack of statistical significance in other scenarios indicates the possible influence of other factors that were not accounted for in this study, such as ocean currents and atmospheric pressure systems.
The critical role of atmospheric and oceanic phenomena in shaping sea ice distribution and variability has been shown before [42][43][44][45][46]. Importantly, the ongoing transition towards thinner, predominantly first-year ice [11,12] could increase the vulnerability of the ASI to atmospheric and oceanic influences, contributing to the observed fluctuations in the sea ice cover.
The growing variability in the ASI area, especially in the context of potentially ice-free summer periods, marks an era of unprecedented environmental change.The fluctuations may manifest as sporadic episodes of extremely low ice coverage, interspersed with phases of temporary ice regrowth, giving the impression of a temporary recovery of the sea ice loss.Moreover, these trends could lead to the onset of ice-free summers sooner than previously anticipated, disrupting environmental and socio-economic balances [3,[26][27][28][29][30]33].Arctic marine mammals, adapted to ice-dependent habitats, would face heightened risks from a seasonally ice-free Arctic, including increased predation, competition, and disease due to northward-expanding temperate species, leading to biodiversity threats and genetic isolation [31,32,34].In parallel, northern Indigenous communities, whose lifestyles and economies depend on stable sea ice conditions, may experience declining food security [33].Additionally, the presence of sea ice and its overlying snow cover plays a crucial role in regulating the amount of sunlight that reaches the upper ocean.This variation in light penetration affects the growth of algae and phytoplankton blooms beneath the ice, further influencing the region's marine ecosystems [47,48].
The increasing fluctuations in the sea ice cover add to the difficulty of predicting future ASI trajectories.To visually represent this concept, we combine a straightforward linear regression prediction approach with the uncertainties derived from our analysis (figures 6 and S19-S21).This is not intended to provide precise future predictions for the ASI, but rather to illustrate how increasing fluctuations contribute to greater uncertainty in projections of the sea ice area.This uncertainty is particularly concerning for navigation.
While the increased open water coverage (tables S8-S11) allows for the establishment of maritime routes across the CA region, the enhanced variability in the sea ice area (figures 4 and S15-S17) may lead to sudden changes in navigability, increasing the risk of delays, rerouting, and potential ice entrapment.As a result, navigators will require more frequent and detailed ice forecasts, increasing operational costs and complexity.
In summary, the anticipated increasing variability in the ASI cover presents numerous challenges.Firstly, it might lead to earlier ice-free summers than previously projected.Secondly, it complicates scientific .To align the model data with the observation-based projection, we shift the variance in time so that the year of complete sea ice loss for the models (2045) aligns with the prediction obtained from the observational data (2064).The red line denotes the threshold for an ice-free Arctic (10 6 km 2 ).forecasting, which in turn impacts navigation, affects geopolitical strategies, and poses challenges for environmental management, calling for a multifaceted approach for adaptation across sectors.

Figure 1 .
Figure 1.ASI area time series and variance averaged over the entire model ensemble under the SSP5-8.5 scenario.(a)The blue solid line is the averaged ASI area from the SSP5-8.5 scenario obtained after synchronizing the model runs to the year of complete ASI loss (see Methods section for details).The filled blue area denotes the values between the 25th percentile (lower line) and the 75th percentile (upper line).The filled grey area corresponds to one standard deviation of the year of complete sea ice loss (2045 ± 11).The red line denotes the threshold for an ice-free Arctic (10 6 km 2 ).(b) same as (a) but for the variance calculated within a window length of w = 50 years.

Figure 3 .
Figure 3. Variance computed for single grid cells under the SSP5-8.5 scenario and averaged over the ensemble.ASI area variance for each grid cell under the SSP5-8.5 scenario calculated within a window length of w = 50 years and normalized at each time step, then averaged over the ensemble (see Methods section for details).(a) 50 years before the complete loss of the ASI, the Central Arctic (CA) sub-region presents low variability.(b) At the moment of complete ASI loss, the highest variance is detected on the northern coast of Greenland and in the area above the Canadian Archipelago, both belonging to the CA sub-region in the MASIE partition (figure 2).

Figure 4 .
Figure 4. Evolution of the probability for ice-free conditions in the CA over 5-year intervals under the SSP5-8.5 scenario.Percentage of the total CA's surface covered by grid cells with different probabilities of being ice-free (defined as sea ice concentration below 15%) over 5-year intervals, averaged over the model ensemble.Fifty years before the total ASI loss, the majority of the CA's surface comprises grid cells having a 0 probability of becoming ice-free in the subsequent five years.After a few decades, we observe a transition to a new, less stable state where most areas are characterized by grid cells having a probability equal to or larger than 0.5 of being ice-free.

Figure 5 .
Figure 5. Difference in surface wind speed variance between models with enhanced ASI variability and models without.(a) Difference in surface wind speed variance computed within a window length of w = 50 years under the SSP5-8.5 scenario 50 years before the complete loss of the ASI (see Methods section for details).The generally positive difference indicated higher surface wind speed variance among the models with increasing ASI area fluctuations, especially over the CA region.(b) Same as (a) but 25 years before the complete loss of the ASI.

Figure 6 .
Figure 6.Increasing uncertainty in ASI area predictions.The gray line represents a linear extrapolation of the observational data for the September ASI area (black dots)[49].The uncertainty is included as the square root of the change in the variance with respect to the years 1961-1990 averaged over the model ensemble under the SSP5-8.5 scenario, calculated within a window length of w = 50 years (figure1(b)).To align the model data with the observation-based projection, we shift the variance in time so that the year of complete sea ice loss for the models (2045) aligns with the prediction obtained from the observational data (2064).The red line denotes the threshold for an ice-free Arctic (10 6 km 2 ).