For a given target unlabeled sample x t , we find the most similar state called the nearest labeled samples  x t ′ , from a library based on distance function. Then, the subsequent evolution x t ′ + τ of the x t ′ are weight averaged based on combination functions to calculate x t +τ . To construct the forecast, the nearest labeled samples are weighted as follows:

where x t +τ is the predicted variable with lead time τ, w i is the weight corresponding to the i ^th selected nearest labeled sample; k is the number of the nearest labeled samples, t ′ is the historical period of t. The w_i values are kept constant in the forecast and do not change with the lead time τ.

A control run named Ice-kNN-Ctrl was constructed according to the traditional kNN model. The traditional kNN model has three main procedures that can influence the predictive skill, namely, the distance function, which measures the similarity between samples, the selection of the k value, and the combination function based on the closest labeled samples (Zhang et al., 2017). Our research focuses on how to organically adapt the physical properties of sea ice to the kNN model, instead of adjusting the parameters. Therefore, Ice-kNN-Ctrl used Euclidean distance to measure the similarity, which is one of the most commonly used distance functions (Zhang et al., 2017); the combination function was set as distance weighting, which assigns weights inverse to the distance and prioritizes the examination of local structures surrounding the samples to be predicted; only SIC was used as input data to calculate the Euclidean distance. A group of hindcast experiments with different k values were conducted and it is found that the prediction results were insensitive to k values ( Figure S1 ). Therefore, the k value for Ice-kNN-Ctrl is only 3.

In this study, following Yang et al. (2020), the summer Arctic sea ice prediction is typically initialized on 1st June. The SIC forecasts are conducted for the independent period from June to September, 2011–2020, which does not overlap with the training period from June to September, 1979 to 2010.

Built on Ice-kNN-Ctrl, we selected different processes for the key steps of the prediction model to optimize its results in Table 1 . These algorithms were identified in advance by our sensitivity experiments to have a considerable impact on the SIC predictive skill. The key steps of the prediction model are: data preprocessing of deseasonalization and detrending; a drift-ice correlation algorithm; expansion of the training library; a distance function; and predictors. Figure 1 illustrates the processes of forecasting the SIC using Ice-kNN.

All experiments were conducted on the Intel Xeon E5-2609 (1.70GHz, 16 cores). The kNN model does not separate the time of training and prediction, so it needs to go through the training and prediction process all over again with each prediction. One 122-day prediction at an initial time costs about 300 seconds.

To assess the forecast skill of the experiments, the SIC predictions of the Ice-kNN model were evaluated at each grid cell using the RMSE of SIC (RMSE_SIC) and bias between the prediction and the observations at 1–122 lead days. The bias is the difference between the prediction and observations for the 10-year average from 2011 to 2020. The Arctic is divided into five regions, as shown in Figure S4 . Owing to the rapid melting of sea ice in recent years, the grid points where SIC have not changed from 2011 to 2020 were excluded when calculating the regional mean predictive skill (Chi and Kim, 2017; Jun Kim et al., 2020).

For SIE verification, three metrics are used: (1) the error of the September monthly mean SIE (Δ^SIE), (2) the RMSE of SIE (RMSE_SIE), and (3) the integrated ice edge error (IIEE) bias Δ^IIEE (Melsom et al., 2019; Fritzner et al., 2020). The total extent (sum of cell areas where SIC > 15%) for each day in September was computed and then averaged for the month for each year into a September mean SIE. The IIEE bias metric, Δ^IIEE, is a measure of the relative difference in sea ice offset predicted by the model. It is computed from three parts, the overestimated and underestimated local SIE and the length of the ice edge. The overestimated part consists of sea ice-free areas that are predicted to be covered with sea ice, and the underestimated part consists of sea ice-covered areas that are predicted to be sea ice free. The length of the ice edge is determined by the ice edge of the observed and predicted fields. A positive Δ^IIEE bias means that the overestimated SIE in the model is large relative to the underestimated SIE, and vice versa.

Sensitivity tests were conducted with the Ice-kNN model using SIC along with one extra variable as a predictor variable (Jun Kim et al., 2020; Liu et al., 2021). To examine the contribution of each predictor to the predictive skill of SIC, the sensitivity is defined as Sens_SIC; to examine the contribution to the predictive skill of SIE, the sensitivity is defined as Sens_SIE. The sensitivity formulas are as follows:

Here RMSE_SIC_sic (RMSE_SIE_sic ) is the RMSE_SIC (RMSE_SIE) of forecast using only SIC, and the RMSE_SIC_predictor (RMSE_SIE_predictor ) is the RMSE_SIC (RMSE_SIE) of forecast using SIC and one extra predictor.

To examine whether the data preprocessing strategy for removing the seasonality and trend could improve the forecast accuracy of the Ice-kNN model, the predictive skill of Ice-kNN-An is compared with that of Ice-kNN-Ctrl in this section. Figure 2 shows the comparison of the hindcast skill between Ice-kNN-An and Ice-kNN-Ctrl measured by the different verification metrics. The RMSE_SIC of Ice-kNN-Ctrl increases with the lead time and stays around 16% after one month ( Figure 2A ). The RMSE_SIC and bias of Ice-kNN-Ctrl in September is mainly distributed in the regions where seasonal sea ice retreats from June to September, including the Beaufort Sea, Chukchi Sea, East Siberian Sea, and Laptev Sea, and it is higher than the climatological prediction in all studied areas ( Figures 3 , 4A ). For SIE, the RMSE_SIE of Ice-kNN-Ctrl increases with lead time and reaches 2.2 × 10⁶ km² in September ( Figure 2B ). Compared with the observations, the prediction of Ice-kNN-Ctrl always overestimates SIE from June to September, and the overestimation gradually increases with the retreat of seasonal sea ice ( Figures 2C , S5A ). It may be related to the difficulty of marginal sea ice prediction and the prediction bias of the Ice-kNN model with the increasing prediction time (Guemas et al., 2016; Liu et al., 2021).

For Ice-kNN-An, there is notable skill enhancement in predicting SIC at lead times longer than one month, and the enhancement of Ice-kNN-An is more pronounced with lead time ( Figure 2A ). According to the spatial pattern, the September RMSE_SIC and bias of Ice-kNN-An in all areas is lower compared with Ice-kNN-Ctrl ( Figures 3 , 4A, B ). For SIE, although both Ice-kNN-An and Ice-kNN-Ctrl tend to overestimate SIE in summer, the predictive skill of Ice-kNN-An is significantly superior than that of Ice-kNN-Ctrl for summer SIE and the improvement increases with lead time ( Figure 2C ). This indicates that the Ice-kNN model can better find the temporal evolution of sea ice by extracting the seasonality and trend, especially from an ice-covered period to an ice-free period.

In extreme ice cover years, such as record low in 2012, the forecast biases are relatively large compared to other years for both models ( Figure 2D ). On the one hand, the lowest minimum Arctic SIE in 2012 is associated with the large multiyear ice volume export and the storm that entered into the central Arctic in early August 2012 (Parkinson & Comiso, 2013; Li et al., 2022). Since the initial day is fixed on June 1st, it is hard for Ice-kNN model to catch the atmospheric disturbance in the extreme cases. On the other hand, for the extreme cases of SIE, Ice-kNN model is not suitable to forecast the extreme values which are not included in the training library due to its prediction principle. While in the other years, Ice-kNN-An shows an impressive improvement compared with Ice-kNN-Ctrl ( Figure 2D ). However, there is still considerable drift ice outside the sea ice edge in both experiments ( Figures S5A, B ).

The improvement of Ice-kNN-An indicates that the deseasonalization and detrending step is useful to improve the Arctic sea ice forecast accuracy of the Ice-kNN model. Therefore, the data preprocessing strategy was used to remove the seasonality and the trend components in subsequent experiments.

Due to the lack of spatial continuity in the traditional kNN model, both Ice-kNN-Ctrl and Ice-kNN-An forecasts showed unrealistic drift ice around the sea ice edge ( Figures S5A, B ). This section studies the impact of different drift-ice correction algorithms in the Ice-kNN prediction. As shown in Figure 5A , the RMSE_SIC of Ice-kNN-F shows poor predictive skill at lead times of less than 30 days, but it is better than Ice-kNN-An at lead times longer than 30 days. According to the distribution of September RMSE_SIC and bias, Ice-kNN-F is superior to Ice-kNN-An in predicting SIC in September in all sea areas except the Chukchi Sea ( Figures 3 , 4B, C ). The results show that the pointwise modeling is better for the prediction of SICA caused by short-term-scale disturbances, whereas for the SICA caused by large-scale anomalies, selecting the similarity using the full-field distance could improve the predictive skill with lead times of more than one month. Ice-kNN-FP predicted the SIC point by point based on Ice-kNN-F determination of the sea ice edge. It performs better than Ice-kNN-F at lead times of less than 30 days and better than Ice-kNN-An at lead times of longer than 30 days ( Figure 5A ). In addition, Ice-kNN-FP has a lower unrealistic drift-ice bias compared with Ice-kNN-An, especially in the Chukchi, East Siberian, Laptev, and Kara seas ( Figures 3 , S5B, D ).

For SIE, the prediction bias of Ice-kNN-F in the short-term lead time is larger. However, Ice-kNN-FP, which combines full-field distance and single-point distance, shows a significant improvement in sea ice edge compared with Ice-kNN-F and Ice-kNN-An for the whole lead time ( Figures 5B, C ). In Figure 5D , it can be seen that Ice-kNN-FP effectively reduces the bias of the September mean SIE in Ice-kNN-F and Ice-kNN-An in most years.

The kNN model, which takes only a single grid point as the prediction sample, lacks physical spatial connection, and leads to the prediction of unrealistic drift ice. A drift-ice correlation algorithm, which selects similarity by full-field distance, would consider the spatial continuity of sea ice but ignore the local SICA caused by short-term disturbance. Therefore, the full-field distance is first used to limit the sea ice coverage, and then pointwise modelling is carried out to predict the SIC of each single grid point, which can effectively correct the unrealistic drift ice of pointwise modelling and the initial sea ice migration bias of full-field modelling. In the following kNN models, the drift-ice correction algorithm of Ice-kNN-FP is applied.

In principle, the nearest neighbors with insufficient similarity could diverge relatively quickly in time compared with the very close nearest neighbors. The limited Ice-kNN forecast skill may therefore be partly due to the relatively small number of available training labeled samples, which makes the nearest neighbor selection a challenge. As the most accurate summer SIC datasets are limited to the satellite era starting in the 1979, the training labeled samples for each state has only 32 to 41 members from 2011 to 2020. To verify the sensitivity of the predictive skill to the number of training labeled samples, Ice-kNN-PA and Ice-kNN-PFA expand the training library by adding adjacent calendar days. The number of training labels samples in Ice-kNN-PA increases from 96 to 123 from 2011 to 2020, and that in Ice-kNN-PFA increases to 123.

Compared with Ice-kNN-FP, the RMSE_SICs of the Ice-kNN-PA and Ice-kNN-PFA seems not to be significantly reduced (Figure 6A), but from the perspective of different sea areas, the Ice-kNN-PA and Ice-kNN-PFA mainly reduced the positive SIC bias in the sea ice marginal zone of Beaufort Sea and the East SiberianLaptev seas (Figures 4D–F, S6A, B). For SIE, Ice-kNN-PA and Ice-kNN-PFA have a larger initial bias at lead times of less than two weeks ( Figure 6B ), which is mainly due to the underestimation of SIE ( Figure 6C ). However, Ice-kNN-PA and Ice-kNN-PFA show significant improvement in predicting SIE after that ( Figures 6B, C ). In Figure 6D , except for 2017 in Ice-kNN-PFA, the September mean SIE bias of Ice-kNN-PA and Ice-kNN-PFA is reduced by about 0.5 million square kilometers compared with Ice-kNN-FP.

In general, expanding the training library will cause an initial underestimation bias, but it will not rapidly diverge with increasing lead time. The predictive skill of the experiments that use the expanded training library are significantly improved compared with Ice-kNN-FP with lead times of more than two weeks. Using the future adjacent calendar days as the training labeled samples has relatively little impact, except for 2017 when Ice-kNN-PFA selects future sample as the nearest sample. The younger and thinner Arctic sea ice in recent years is more sensitive to external forcing (Parkinson & Comiso, 2013), resulting in a large deviation in the forecast results when the future adjacent calendar days are selected in the training library.

Therefore, selecting sufficient training labeled samples could better improve the predictive skill of the Ice-kNN model for Arctic sea ice. In subsequent experiments, the strategy of expanding the training library by adding past adjacent calendar days as training labeled samples to predict the Arctic sea ice was applied.

In previous Ice-kNN models, the Euclidean distance has been selected most frequently as the distance function. In the prediction of sea ice, not only the distance between grid cells but also the spatial correlation coefficient between states should be considered to select similarity. In this section, a compound distance function scheme, including the spatial anomaly correlation coefficient and the Euclidean distance between sea ice, is studied.

As shown in Figure 7A , the RMSE_SIC of Ice-kNN-PC is comparable to that of Ice-kNN-PA. According to the spatial pattern, the September RMSE_SIC and bias of Ice-kNN-PC decreases in the East Siberian–Laptev seas and Kara–Barents–Greenland seas, but increases in the Beaufort Sea, Chukchi Sea, and Baffin Bay–Canadian Archipelago compared with Ice-kNN-PA ( Figures 3 , 4F, G ). For the SIE, the RMSE_SIE of Ice-kNN-PC has a larger bias at lead times of less than one week ( Figure 7B ), which is mainly due to the underestimation of SIE ( Figure 7C ). However, Ice-kNN-PC shows improvement in SIE compared with Ice-kNN-PA at lead times of more than one week and the improvement is more pronounced with increasing lead time ( Figures 7B, C ). In the Figure 7D , except for 2012, the biases of monthly mean SIE in September of Ice-kNN-PC from 2011 to 2020 are within about 0.5 million square kilometers, which is lower than Ice-kNN-PA.

In general, the composite distance function with a spatial anomaly correlation coefficient is beneficial to the prediction of Arctic sea ice. The new distance function considers not only the similarity between samples at a single point through the Euclidean distance but also the spatial mode of the SICA through the spatial anomaly correlation coefficient. It is helpful for the Ice-kNN model to consider the large-scale spatial variation of SICA when selecting the similarity. Therefore, in the following experiments, a composite distance function combining the Euclidean distance and the spatial anomaly correlation coefficient is used to further improve the predictive skill of the Ice-kNN model.

To verify the impact of sea ice-related atmospheric and oceanic variables on the predictive skill of summer Arctic sea ice, SIC and SIE sensitivity indices including sea ice-related variables were calculated based on the Ice-kNN-PC model in Figure 8 . A positive sensitivity index indicates that a variable has a positive contribution to the predictive skill of SIC (SIE).

Ice-kNN-Sflux, which selects Sflux and SIC as predictors, improves the predictive skill of SIC for the whole lead time ( Figure 8A ). SST improves the predictive skill at lead times of less than one month, but SLP and T2m provide only limited improvement in the predictive skill at lead times of about one to two months ( Figure 8A ). According to the distribution of RMSE_SIC and bias in September, the improvement of the predictive skill of Ice-kNN-Sflux mainly occurs in the Beaufort Sea, compared with Ice-kNN-PC ( Figures 3A , 4G–K ). For SIE, the sensitivity index is calculated based on the RMSE_SIE. All the sea ice-related variables show improvement of the SIE predictive skill at lead times of less than one month. But for lead times longer than one month, all the sea ice-related variables show negative contributions to the predictive skill of SIE.

The predictive skill of Ice-kNN-Sflux for SIC is significantly better than both the climatological prediction and the anomaly persistence prediction at lead times of longer than two weeks ( Figure S7 ). The predictive skill of Ice-kNN-Sflux for September SIC is significantly better than the anomaly persistence prediction for the whole Arctic, and significantly better than the climatological prediction for the whole Arctic, except for Baffin Bay and the Canadian Islands ( Figure 3 ). The prediction bias of Ice-kNN-Sflux in September SIE is reduced by 2.0 × 10⁶ km² compared with the climatological prediction and by 3.0 × 10⁶ km² compared with the anomaly persistence prediction.

In general, for daily Arctic sea ice forecasts in summer, the Sflux fields, which have a direct relation to sea ice (Liu et al., 2021), can enhance the predictive skill of sea ice. SLP and T2m show little improvement of the predictive skill of sea ice, which may result from the chaotic behavior of the atmosphere (Mohammadi-Aragh et al., 2018). While the surface oceanic field SST does not show its long-term memory, which may result from the interpolation bias. Similar results were obtained in Liu et al. (2021) using the deep learning model ConvLSTM.

The Sea Ice Prediction Network (SIPN) is an open platform that has been collecting predictions of Arctic SIE in September around the world since 2008, then compiling and presenting them to those interested in Arctic sea ice. September SIE predictions have been submitted to the SIPN in June, July, and August since 2008, with an additional September submission added in 2021. There are a variety of prediction methods, including heuristic, statistical, mixed, dynamic, and machine learning/other. To further evaluate the Arctic sea ice forecast skill of Ice-kNN model, we compared the Arctic September SIEs in 2021 and 2022 using the Ice-kNN-Sflux model with the observations and the contributions for the September SIE predictions to the Post-Season Sea Ice Outlook for 2021 and 2022 (Bhatt et al., 2022a; Bhatt et al., 2022b). We utilized the SIC of May 31, June 30, July 31 and August 31 respectively as the input of Ice-kNN-Sflux model to predict the September SIEs. It should be noted here we use a hindcast (not real-time forecast) result of the Ice-kNN-Sflux model.

As shown in Figure 9A , the observed September SIE in 2021 was 4.92 × 10⁶ km² (reported by NSIDC). The median hindcasting result of Ice-kNN-Sflux from June to September is 4.8 × 10⁶ km², with a quartile range of 4.62 to 5.04 × 10⁶ km² ( Figure 9A ). The median September estimate based on all contributors of SIPN were 4.37 × 10⁶, 4.36 × 10⁶, 4.39 × 10⁶, and 4.39 × 10⁶ km², respectively, from June to September. In comparison, our hindcasts using Ice-kNN-Sflux were 5.34 × 10⁶, 4.94 × 10⁶, 4.66 × 10⁶, and 4.49 × 10⁶ km², respectively ( Table S1A ). For 2022, the medians September estimate of SIPN were 4.57 × 10⁶, 4.64 × 10⁶, 4.83 × 10⁶, and 4.91 × 10⁶ km², respectively, from June to September, which approaches observation 4.87 × 10⁶ km² (reported by the NSIDC). The hindcasts of Ice-kNN-Sflux from June to September were 5.05 × 10⁶, 4.47 × 10⁶, 5.65 × 10⁶, and 4.62 × 10⁶ km² ( Table S1B ). The median hindcasting result of Ice-kNN-Sflux from June to September is 4.835 × 10⁶ km², with a quartile range of 4.58 to 5.2 × 10⁶ km² ( Figure 9B ).