Edge displacement scores

. As a consequence of a diminishing sea ice cover in the Arctic, activity is on the rise. The position of the sea ice edge, which is generally taken to deﬁne the extent of the ice cover, changes in response to dynamic and thermodynamic processes. Forecasts for sea ice expansion due to an advancing ice edge will provide information that can be of signiﬁcance for operations in polar regions. However, the value of this information depends on the quality of the forecasts. Here, we present methods for examining the quality of forecasted sea ice expansion and the geographic location where the largest expansion are 5 expected from the forecast results. The algorithm is simple to implement, and an examination of two years of model results and accompanying observations demonstrates the usefulness of the analysis.


Introduction
Due to climate change the sea ice extent is in decline in the Arctic (Parkinson, 2014). This change has led to increased activity 10 in the region, and commercial shipping in open waters via Arctic sea routes will become increasingly economically viable in the 21 st century (Aksenov et al., 2017). Thus, products for monitoring and forecasting sea ice conditions are receiving growing attention.
The past years have seen a flurry of activity related to assessing the quality of sea ice products. Dukhovskoy et al. (2015) present a review and comparison of various traditional metrics for assessments of the skill of sea ice models. Goessling et 15 al. (2016) introduce the Integrated Ice-Edge Error (IIEE), a quantity for describing mismatching sea ice extents from two products, in their examination of the predictability of the sea ice edge. Melsom et al. (2019) took advantage of the IIEE in their examination of various metrics for assessment of the quality of forecasts for the sea ice edge position. Methods for examining the quality of probabilistic results for sea ice conditions have been introduced by Goessling and Jung (2018) and . Recently, Cheng et al. (2020) have examined the accuracy of a visually estimated ice concentrations monitoring 20 product.
The changing position of the sea ice edge is generally not only shifted by dynamic advection, but can be significantly affected by the thermodynamics as well (Bitz et al., 2005). Thus, the temporal displacement of the sea ice edge will be affected by freezing along the perimeter of the sea ice extent in winter, and melting in summer. Hence, pattern-recognition algorithms for displacements using maximum cross-correlation (MCC) methods such as those introduced by Leese et al. (1971) for wind vectors, and later for ocean surface currents (Tokmakian et al., 1990) and sea ice vectors (Lavergne et al., 2010), are not ideal for tracking displacements of the sea ice edge. Ebert and McBride (2000) examined the position error of the contiguous rain area in weather forecasts. They determined the error vector from a total squared error minimization method when shifting the forecasted rain region to match the corresponding observations. Their preference of applying an error minimization algorithm rather than an MCC approach was motivated by 30 the former having better representations of displacement of rainfall maxima. Displacement of the perimeter of the contiguous rain area was not addressed in their investigation.
We begin this study by presenting a new algorithm for assessing the quality of representations of the sea ice edge by comparing results for two different products. This is described in Sect. 2, where issues related to sea ice emerging in the vicinity of boundaries are also addressed. Next, in Sect. 3 we apply the algorithm in an examination of displacements of the 35 sea ice edge in the Barents Sea in a model product, and compare the results to data from an observational product. Finally, we provide our concluding remarks in Sect. 4.

Methods
In order to illustrate the validation metrics that are introduced in this Sect., some idealized distributions are introduced, as depicted in Fig. 1. The domain is divided into 1000 × 500 square grid cells, and we set the length of the side of a grid cell 40 to 1. Denote the line that separates regions with binary values 0 and 1 as an edge line, and let L (o) (t) and L (m) (t) denote observed and modeled edges, respectively, at time t. Idealized examples with edges for L (o) and L (m) at two different times, t 0 and t 0 + ∆t, are displayed. In the context of forecasting, L (m) (t 0 ) may be taken to represent the model initialization at t 0 and L (m) (t 0 +∆t) is then the forecast at a temporal range of ∆t. The other binary fields can represent observations at the same times.

Single product metrics
We aim at defining metrics that describe differences in maximum edge displacements between two products. In order to do so, we must first introduce a quantity that properly measures the maximum displacement in one product. Here a definition is provided which is a gridded, signed, one-sided variation of the Hausdorff distance (Dukhovskoy et al., 2015).
For the remainder of this investigation we will take the binary fields to be representations of sea ice, with values assigned 50 to 0 and 1 for conditions of no ice and ice, respectively. We will here associate the presence of ice (value 1) with sea ice concentration c exceeding c edge = 0.15. In a gridded representation the ice edge can then be taken to be constituted by the grid cells e = [i, j] that meet the condition where ∧ is the logical AND operator. Denoting the N grid cells that satisfy this condition by e 1 , e 2 , . . . , e N the ice edge is then 55 the line This follows the algorithm presented in Melsom et al. (2019). Let L (1) , L (2) denote the sea ice edges for two representations of the sea ice cover. Furthermore, let d 2:1 n be the displacement distance between grid node e (2) n in L (2) and line L (1) , i.e., where s is +1 or -1 when e (2) n is on the no ice or ice side of L (1) , respectively, i.e., s is +1 if c[e  (1) is the sea ice concentration for grid cell e for the representation with an ice edge given by L (1) ). Here, ||z|| is the Euclidean distance of z. We can now introduce the maximum distance as so that d 2:1 max will represent the displacement of the largest sea ice advance from L (1) to L (2) . For reference, we note that the Hausdorff distance d H between lines L (2) and L (1) is see e.g. Dukhovskoy et al. (2015).

70
So far, we have restricted the analysis to consider the maximum displacement. However, it is of interest to examine more generally the displacement distance between L (2) and L (1) . In addition to a visual inspection that can be done by looking at a graphical presentation (see Fig. 1 for examples), we can consider all values for d 2:1 n and summarize the distribution in a table with a suitable definition of distance categories, or present the distribution as e.g. a cumulative probability distribution.
The distribution of selected distance categories for the idealized model results displayed in Fig. 1   To avoid inflating the sample size when time series results are examined, one may consider subsampling at the spatial decorrelation length along the ice edge. For the distribution of d 2:1 n a proper decorrelation length can be computed if the edge nodes e (2) n are in sequence along L (2) .

80
Also, when a time series of results is examined, the distribution of d 2:1 max (t) can be examined analogously to results for d 2:1 n . We will show results for d 2:1 max (t) when comparing model results and observations in the case study that is undertaken in Sect. 3.

Two-product metrics
The main purpose of the work presented here is to define metrics that can contribute in an evaluation of the quality of model  From the perspective of an observer, a useful property is the quality of the forecasted maximum displacement of the binary field, over the forecast period. A simple quantity that provides relevant information, is the difference in the maximum displacement as given by Eq. (4), i.e., where o1, o2 are observed ice edges at t 0 and t 0 + ∆t, respectively (black and gray lines in Fig. 1), and m1, m2 are the corresponding model results. From the results for the idealized example that was introduced in Sect. 2.1 above the model is over-estimating the maximum displacement, by ∆d 2:1 max = 15.3 grid cell units. A similar quantity that provides site specific information is the local difference in displacement of the model binary field where the maximum value is found in the observations. Let e o2 0 be the position in L (o2) to which the maximum edge displace- where L (m1) = L (m) (t 0 ). For the idealized example, we find that δ m2:m1 0 = 83.9. The local difference in displacement between 105 model and observations, with reference to the position e o2 0 , becomes so, for the idealized example we have ∆δ 2:1 max = −14, i.e. a local underestimation of the displacement in the model results. One aspect which is not disclosed by the metrics introduced thus far, is to what degree forecasts manage to reproduce the geographical location of the observed maximum displacements. In order examine such a relation, we first compute the 110 decorrelation length of displacements given by Eq. (3). If we denote this grid distance by ∆n, we restrict the analysis of grid cells and corresponding displacements to {. . . , δ m2:m1 0−2∆n , δ m2:m1 0−∆n , δ m2:m1 0+∆n , δ m2:m1 0+2∆n , . . .} respectively, limited by the first and last nodes along the line L (m2) . Next, we construct bins analogously to the method used for 115 producing rank histograms (Talagrand diagrams) for ensemble forecasts (Hamill, 2001): First, distances listed in Eq. (10) are sorted by increasing values, and then bins are introduced for values smaller than the minimum distance, the intervals between the sorted distances, and for values larger than the maximum distance. The bin placement of δ m2:m1 0 then gives the rank of this displacement.
In the present idealized example we find that ∆n = 42, and the rank of δ m2:m1 0 in the 24 resulting bins is 9. When multiple 120 forecasts are examined, the decorrelation length will generally change, as will the length of the edges. Thus, in order to derive a meaningful statistic quantity we subsample a fixed sized random set of grid cells from Eq. (10), and an analysis of the ranks of displacement distances can be performed.
For the idealized example, a set of nine randomly subsampled edge positions from those given by Eq. (9) for the model results at t = t 0 + ∆t is displayed by open circles in Fig. 1. For this particular case, in the range from 1 to 10 the rank of the 125 displacement δ m2:m1 0 is 3.

Open boundaries and coasts
Sea ice products may be regional, where e n OB is any ocean grid cell along the boundary of the domain. Then L (1) can be replaced bỹ 135 and for the corresponding distances we introduce the notationd, so e.g. Eq.
(3) is writteñ It must be noted that if the ice is imported into the domain, the distancesd associated with such a displacement will be underestimated, since the real position of the ice edge outside of the analysis domain at t 0 is unknown.
Similarly, there can be cases where freezing of ice occurs along the coastline, e.g. due to colder air in the vicinity of 140 continents, or less salty water masses close to the coastline. This is another case where the algorithm above can yield grossly exaggerated displacement distances. Again, the problem can be overcome by including additional grid lines.
Set the coastal grid lines as where e n C is any ocean grid cell along the coastline. Then L (1) can be replaced by For a regional model, the typical situation is that there are both open boundaries and coastlines. In that case, we may combine the above modifications of the algorithm by adopting  (Lien et al., 2013).
The ocean module of the coupled model used for the regional simulation is the Regional Ocean Modeling System (ROMS), 155 described in Haidvogel et al. (2008) and references therein. The sea ice module was developed by Budgell (2005). The ice model dynamics are based on the elastic-viscous-plastic (EVP) rheology after Hunke and Dukowicz (1997) and Hunke (2001), and the ice thermodynamics are based on Mellor and Kantha (1989) and Hakkinen and Mellor (1992). The model results for sea ice concentration are somewhat noisy on the grid cell scale, owing to the dispersiveness of the numerical scheme. In some regions, the grid cells that constitute the ice edge as defined by Eq. (1) can then appear as a mesh-160 like collection of cells. In order to reduce the impact of this issue, we applied the second order checkerboard suppression algorithm (Li et al., 2001) to the model results before conducting the present analysis.
We compare model results with observations from the Arctic Ocean Sea Ice Concentration Charts Svalbard which is a multi-sensor product that uses data from Synthetic Aperture Radar (SAR) instruments as its primary source of information (WMO, 2017). This product covers the northern Nordic Seas, the Barents Sea and adjacent ocean regions. It is available on 165 a stereographic grid projection with a resolution of 1 km. The product will be referred to as the ice chart data hereafter. Data availability is restricted to working days. During a regular week, we then have four days with 24 h displacement results. The data set is also slightly reduced due to holidays, and a total of 354 days with 24 h ice edge displacement results were available from the present two year period.
The present study will be restricted to results and data for the Barents Sea. The SVIM simulation domain is displayed in   We first examine the distribution of daily maximum ice edge displacements. The results are summarized in Table 2, where results from the 354 days with 24 h displacements from both products have been included. We note that about 2/3 of the 175 displacements in model results are in the range 10 -30 km. The corresponding distribution of results from the ice chart data has two maxima, one for the range 20 -40 km which accounts for nearly half of the cases, and a secondary maximum for short (0 -10 km) displacements. The averages of the daily maximum displacement distances are 25 km and 36 km for the SVIM results and the ice chart data, respectively.
The category distributions in Table 2  The resulting frequency distribution for each of the ten ranks is displayed as gray vertical bars in Fig. 5, with rank values from 0 to 9. The highest rank (9) results when the model displacement close to the site with maximum displacement in the observations (the reference displacement, δ m2:m1 0 ) is larger than all displacements from the nine subsampled ice edge positions.

195
The next rank (8) corresponds to cases where one and only one of the subsampled positions have a larger displacement than the reference displacement, and so on. In other words, high ranks indicate situations in which the position of the maximum displacement is described with a relatively high quality.
The expected average rank from a model with no quality in reproducing the position of the maximum ice edge displacement is 4.5, with half of the cases with a rank in the range 0 -4, and the other half in the range 5 -9. In the analysis of the results 200 from the SVIM archive, we find that 35% of the cases have rankings in the range 0 -4, while 65% have rankings in the range 5 -9.
Furthermore, the average rank in the present analysis is 5.5. For a random distribution of 235 integer numbers in the range 0 -9 the 0.005 th and 0.995 th percentiles are 4.015 and 4.985, respectively. Thus, the analysis reveals that while the model results are far from perfect, the average ranking of 5.5 is significantly higher than results from random spatial distributions of ice edge 205 displacements.

Concluding remarks
We present a simple algorithm for examination of the displacement of the edge (or the front) of a binary field. This forms the basis for a subsequent analysis of statistical properties for such displacements. Furthermore, additional methods have been introduced for the purpose of comparing two different products that are available as descriptions of the same situation.

210
The present study has been framed in the context of results for displacements of the sea ice edge. Thus, the case study which was investigated in Sect. 3 was based on data for the sea ice edge from satellite observations, and simulation results from a coupled ocean -sea ice model. However, the algorithm that was introduced in Sect. 2 can be applied to the displacement of the perimeter of any property that can be represented by a continuous binary field. Stratiform precipitation is an example of another property for which the present methods are relevant.

215
Note that we have used the term displacement rather than advection. The reason for this is that displacements need not be purely of an advective nature. In the case of sea ice, the displacement of the initial edge will generally be affected by freezing or melting along the perimeter of the sea ice extent. Analogously, displacement of the area affected by stratiform precipitation can be affected by new condensation or partial depletion of the cloud.
As demonstrated in the example depicted in Fig. 4 as revealed from the sample situation in Fig. 4.
However, there may be other issues that can distort results that are produced by the present analysis. One example is cases 225 where features are seen to arise seemingly spontaneous from one time of analysis to another: The algorithm in Sect. 2 can e.g., if applied to precipitation data, give rise to unrealistic results for displacements when convective precipitation cells become established.
Results from the algorithms that are introduced in the present study give valuable information regarding the changing extent of sea ice, and how well the displacement of the observed and modeled sea ice edge agree. These algorithms have proven to 230 provide simple, yet robust and informative assessments for the quality of ice edge forecasts both with respect to the largest displacements from one time to another as well as with respect to the reproduction of the geographical position where the largest displacement occurs.
Code and data availability. The idealized distributions of concentrations that are depicted in Fig. 1, and source code for computing results displayed in Table 1 and Fig. 2, are available as https://doi.org/10.5281/zenodo.4545686 (Melsom , 2021