The size distribution of ice floes in the polar seas affects the dynamics and thermodynamics of the ice cover, and ice-ocean models are beginning to include the floe size distribution (FSD) in their simulations. The FSD has previously been reported to follow a power law of the form *x ^{–α}*, where

*x*is the floe size and –

*α*characterizes how steeply the number of floes decreases as

*x*increases. Different studies have found different values of

*α*and different ranges of

*x*over which the power law applies. We found that a power law describes the FSD in the Beaufort and Chukchi seas reasonably well over floe sizes from 2 to 30 km, based on 187 visible-band satellite images (resolution 250 m) acquired during spring through fall of 2013 and 2014. The mean power-law exponent goes through a seasonal cycle in which

*α*increases in spring, peaks in July or August, and decreases in fall. June is the transition month from spring FSD to summer FSD. This cycle is consistent with the processes of floe break-up in spring followed by preferential melting of smaller floes in summer and the return of larger floes after fall freeze-up. We also analyzed 12 high-resolution satellite images acquired near the low-resolution images in space and time. We found that the FSDs from the high-resolution images follow power laws over floe sizes from 10 m to 3 km. While the power-law exponents of the corresponding high- and low-resolution images do not always match in a strict statistical sense, they suggest the plausibility that the FSD follows a single power law over a wide range of floe sizes. This study covers a larger spatial and temporal sampling space and is based on more satellite images than previous studies. Results have been used for model calibration and validation.

## 1. Introduction

Large parts of the Arctic sea-ice cover consist of distinct ice floes, either frozen together or separated by open water. The sea-ice floe size distribution (FSD) is the number of floes in different size categories in a region, divided by the area of the region (Rothrock and Thorndike, 1984). This is also known as the floe number density. An alternative description of the FSD is the fractional area covered by floes in different size categories, or the floe area density. The FSD is important for several reasons. As small floes melt faster than big floes due to their larger surface area per unit volume, knowledge of the FSD is needed for calculating total ice melt (Steele, 1992). The FSD also affects sea-ice dynamics through the transmission of internal stresses in the ice pack (Zhang et al., 2015). Sea-ice roughness and the transfer of energy and heat between the atmosphere and ocean are affected by the FSD and the presence of floe edges (Birnbaum and Lupkes, 2002; Tsamados et al., 2014; Martin et al., 2016). Thus the partition of floes into size categories is needed to represent a variety of processes more accurately in models, some of which are beginning to implement FSDs (Horvat and Tziperman, 2015; Hosekova et al., 2015; Zhang et al., 2015, 2016). These models use the floe area density while most observational studies use the floe number density, but the two are mathematically related and there is no complication in switching from one to the other. In this work we consider only the floe number density.

The goal of this work is to characterize the FSD in the Beaufort and Chukchi seas; i.e., to find the distribution of floe sizes, and how that distribution changes over the course of the spring–summer–fall. We are motivated by a recent Office of Naval Research initiative to study the processes that govern the evolution of the Marginal Ice Zone (MIZ) over the course of the year, from spring break-up through fall freeze-up (Lee et al., 2012). The MIZ is the transition region from open water to pack ice, characterized by sea-ice floes of varying shapes and sizes. The present study applies to the MIZ as well as those parts of the pack ice that consist of distinct floes.

Previous studies (Table 1) have reported that the FSD follows a power-law relationship of the form *n*(*x*) = *cx ^{–α}*, where

*x*is the size of a floe and

*n*(

*x*) is the number of floes of size

*x*per unit area. The exponent –

*α*characterizes the shape of the distribution;

*c*is a normalization constant. In a satellite image or aerial photo, the range of

*x*for which the power-law relationship is valid depends on the spatial resolution of the image (at the small end of the size range) and the size of the image itself (at the large end of the size range). The distinguishing feature of a power-law size distribution is its scaling property – there is no natural length scale, so “zooming in” produces an image that looks statistically like the original (Mandelbrot, 1982).

Study . | Location . | Time period . | Cum. or Non-cum.^{a}
. | Exponent^{b}
. | Floe size range (m) . | Number of images^{c}
. |
---|---|---|---|---|---|---|

Rothrock and Thorndike (1984) | Beaufort Sea | June–Aug 1973–1975 | Cum. | –1.7 to –2.5 | 100 to 20,000 | 7 |

Matsushita (1985) | Sea of Okhotsk | not given | Cum. | –2.2 | not given | 1 |

Kergomard (1989) | Fram Strait | June 1987 | Cum. | –1.0 to –1.8 as ice edge is approached | 5 to 2000 | 5 |

Lensu (1990) | Weddell Sea, Antarctica | Feb 1990 | Cum. | –1.4 | 1 to 100 | 1 |

Holt and Martin (2001) | Beaufort and Chukchi seas | Aug 1992 | Cum. | –1.9 to –2.2 (Bfort) –2.2 to –2.6 (Chuk) | 900 to 10,000 | 94 |

Paget et al. (2001) | East Antarctica | Aug 1995 | Non-Cum. | –1.9 to –3.5 as ice edge is approached | 1 to 150 | 6 |

Inoue et al. (2004) | Sea of Okhotsk | Feb 18, 2000 | Cum. | –1.5 to –2.1 as ice edge is approached | 7 to 53 | 2 |

Toyota et al. (2006) | Sea of Okhotsk | Feb 2003 | Cum. | –1.2x < 40 m–1.9 x > 40 m | 1 to 1000 | 3 |

Toyota et al. (2011) | Weddell Sea, Antarctica | Sept–Oct 2006, 2007 | Cum. | –1.0 to –1.5x < 40 m–3.2 to –7.6 x > 40 m | 2 to 100 | 122 |

Toyota et al. (2016) | Wilkes Land, East Antarctica | Sept–Nov 2012 | Cum. | –1.3 to –1.4x < 100 m–2.9 to –3.1 x > 1500 m | 5 to 100 and 1500 to 10,000 | 16 |

Steer et al. (2008) | Weddell Sea, Antarctica | Dec 2004 | Non-cum. | –1.9x < 20 m–2.8 to –3.4 x > 20 m | 2 to 400 | 130 |

Lu et al. (2008) | Prydz Bay, East Antarctica | Dec 2004 to Feb 2005 | Cum. | –0.6 to –1.4 as ice edge is approached | 2 to 100 | 19 |

Perovich and Jones (2014) | Beaufort Sea | June–Sept 1998 | Cum. | –2.0 to –2.2 following a seasonal cycle | 10 to 10,000 | 10 |

Gherardi and Lagomarsino (2015) | Barents/Kara seas and Weddell Sea | Spring 2000, 2001, 2003, 2009 | Non-cum. | –2.0 | 2 to 5000 | 4 |

Geise et al. (2016)^{d} | East Siberian Sea | June–Aug 2000–2002 | Cum. | –0.6 to –1.2x < 200–850 m–1.2 to –2.0 x > 200–850 m | 12 to 3900^{d} | 6 |

Wang et al. (2016) | Beaufort and Chukchi seas | Summer–fall 2014 | Cum. | –1.0 to –1.5 | 5 to 10,000 | 18 |

Hwang et al. (2017a) | Beaufort Sea | Aug–Sept 2014 | Non-cum. | –2.7 to –3.0 | 450 to 3000 | 29 |

This study | Beaufort and Chukchi seas | Mar–Oct 2013, 2014 | Non-cum. | –1.9 to –2.8 following a seasonal cycle | 10 to 30,000 | 273 |

Study . | Location . | Time period . | Cum. or Non-cum.^{a}
. | Exponent^{b}
. | Floe size range (m) . | Number of images^{c}
. |
---|---|---|---|---|---|---|

Rothrock and Thorndike (1984) | Beaufort Sea | June–Aug 1973–1975 | Cum. | –1.7 to –2.5 | 100 to 20,000 | 7 |

Matsushita (1985) | Sea of Okhotsk | not given | Cum. | –2.2 | not given | 1 |

Kergomard (1989) | Fram Strait | June 1987 | Cum. | –1.0 to –1.8 as ice edge is approached | 5 to 2000 | 5 |

Lensu (1990) | Weddell Sea, Antarctica | Feb 1990 | Cum. | –1.4 | 1 to 100 | 1 |

Holt and Martin (2001) | Beaufort and Chukchi seas | Aug 1992 | Cum. | –1.9 to –2.2 (Bfort) –2.2 to –2.6 (Chuk) | 900 to 10,000 | 94 |

Paget et al. (2001) | East Antarctica | Aug 1995 | Non-Cum. | –1.9 to –3.5 as ice edge is approached | 1 to 150 | 6 |

Inoue et al. (2004) | Sea of Okhotsk | Feb 18, 2000 | Cum. | –1.5 to –2.1 as ice edge is approached | 7 to 53 | 2 |

Toyota et al. (2006) | Sea of Okhotsk | Feb 2003 | Cum. | –1.2x < 40 m–1.9 x > 40 m | 1 to 1000 | 3 |

Toyota et al. (2011) | Weddell Sea, Antarctica | Sept–Oct 2006, 2007 | Cum. | –1.0 to –1.5x < 40 m–3.2 to –7.6 x > 40 m | 2 to 100 | 122 |

Toyota et al. (2016) | Wilkes Land, East Antarctica | Sept–Nov 2012 | Cum. | –1.3 to –1.4x < 100 m–2.9 to –3.1 x > 1500 m | 5 to 100 and 1500 to 10,000 | 16 |

Steer et al. (2008) | Weddell Sea, Antarctica | Dec 2004 | Non-cum. | –1.9x < 20 m–2.8 to –3.4 x > 20 m | 2 to 400 | 130 |

Lu et al. (2008) | Prydz Bay, East Antarctica | Dec 2004 to Feb 2005 | Cum. | –0.6 to –1.4 as ice edge is approached | 2 to 100 | 19 |

Perovich and Jones (2014) | Beaufort Sea | June–Sept 1998 | Cum. | –2.0 to –2.2 following a seasonal cycle | 10 to 10,000 | 10 |

Gherardi and Lagomarsino (2015) | Barents/Kara seas and Weddell Sea | Spring 2000, 2001, 2003, 2009 | Non-cum. | –2.0 | 2 to 5000 | 4 |

Geise et al. (2016)^{d} | East Siberian Sea | June–Aug 2000–2002 | Cum. | –0.6 to –1.2x < 200–850 m–1.2 to –2.0 x > 200–850 m | 12 to 3900^{d} | 6 |

Wang et al. (2016) | Beaufort and Chukchi seas | Summer–fall 2014 | Cum. | –1.0 to –1.5 | 5 to 10,000 | 18 |

Hwang et al. (2017a) | Beaufort Sea | Aug–Sept 2014 | Non-cum. | –2.7 to –3.0 | 450 to 3000 | 29 |

This study | Beaufort and Chukchi seas | Mar–Oct 2013, 2014 | Non-cum. | –1.9 to –2.8 following a seasonal cycle | 10 to 30,000 | 273 |

^{a} Cumulative (Cum.) or non-cumulative (Non-cum.) representation of FSD.

^{b} Exponent (or slope in log-log space) of power law describing FSD.

^{c} Number of images, photographs, or frames used to derive FSD. When multiple frames were stitched together to form a larger composite image from which one FSD was derived, only one image was counted. When several regions from the same image were analyzed separately, each region was counted as a separate image.

^{d} Geise et al. (2016) reported cumulative FSD vs. floe area (*A*), not floe diameter (*x*). As *A* = *kx*^{2}, where *k* is a constant, we have multiplied their reported exponents by 2 for this table. Their reported floe areas range from 10^{2} to 10^{7} m^{2}. We have converted these to floe sizes 12 to 3900 m using *k* = 0.66 (from Rothrock and Thorndike, 1984). Their reported transition scale separating two power-law regimes is 28 × 10^{3} to 485 × 10^{3} m^{2}. We have converted these to 200–850 m for this table (again using *k* = 0.66).

In this work we show that the FSD in the Beaufort and Chukchi seas follows a power law whose exponent undergoes a seasonal transformation during spring and summer that is consistent with floe break-up and melting. A companion paper to this one, *On reconciling disparate studies of the sea-ice floe size distribution* (Stern et al., 2018), addresses the representation of FSDs and compares results from different investigations.

## 2. Data

The Moderate Resolution Imaging Spectroradiometer (MODIS) is a 36-band instrument on NASA’s Terra and Aqua satellites (http://modis.gsfc.nasa.gov/), providing daily coverage of the Beaufort and Chukchi seas. We used the visible band images (band 3), which have a pixel size of 250 m. We downloaded the images automatically from the NASA Global Imagery Browse Service via the Web Map Tile Service (WMTS) interface (http://map1.vis.earthdata.nasa.gov/wmts-arctic/MODIS_Terra_CorrectedReflectance_TrueColor/). The images are geolocated (in geoTiff format) and are about 2000 km on a side. Because of clouds, some images are not usable. For images with substantial cloud-free regions, we extracted and analyzed those regions (see Methods). We have analyzed 116 cloud-free regions from 2013 and 140 regions from 2014 (Table 2).

Year . | Number of images . | Number of regions . | Area (10^{3} km^{2})
. | Number of floes . | ||||
---|---|---|---|---|---|---|---|---|

Min . | Median . | Max . | Min . | Median . | Max . | |||

2013 | 83 | 116 | 2.4 | 34 | 360 | 128 | 2563 | 33,716 |

2014 | 104 | 140 | 2.9 | 28 | 130 | 179 | 1587 | 10,586 |

Year . | Number of images . | Number of regions . | Area (10^{3} km^{2})
. | Number of floes . | ||||
---|---|---|---|---|---|---|---|---|

Min . | Median . | Max . | Min . | Median . | Max . | |||

2013 | 83 | 116 | 2.4 | 34 | 360 | 128 | 2563 | 33,716 |

2014 | 104 | 140 | 2.9 | 28 | 130 | 179 | 1587 | 10,586 |

Very high resolution visible images of Arctic sea ice in the Beaufort Sea are available from the USGS Global Fiducials Library (http://gfl.usgs.gov/) as a result of the declassification effort of the MEDEA group (Kwok and Untersteiner, 2011) (hence we refer to them as MEDEA images). The images have a pixel size of 1 m and range from about 10 to 50 km on a side. We have analyzed subsets of three MEDEA images in detail that overlap with MODIS images and were acquired within 3 days of them (Table 3).

Date . | Quantity . | Subset #1 . | Subset #2 . | Subset #3 . | Subset #4 . | Subset #5 . |
---|---|---|---|---|---|---|

May 30 | Area (km^{2}) | 64 | ||||

Num. of floes | 3643 | |||||

Slope of FSD | –2.23 | |||||

June 17 | Area (km^{2}) | 101 | 95 | |||

Num. of floes | 4211 | 2928 | ||||

Slope of FSD | –2.15 | –2.09 | ||||

July 8 | Area (km^{2}) | 16 | 41 | 48 | 52 | 26 |

Num. of floes | 1460 | 4075 | 2322 | 2169 | 1011 | |

Slope of FSD | –2.13 | –2.13 | –2.17 | –2.15 | –2.10 |

Date . | Quantity . | Subset #1 . | Subset #2 . | Subset #3 . | Subset #4 . | Subset #5 . |
---|---|---|---|---|---|---|

May 30 | Area (km^{2}) | 64 | ||||

Num. of floes | 3643 | |||||

Slope of FSD | –2.23 | |||||

June 17 | Area (km^{2}) | 101 | 95 | |||

Num. of floes | 4211 | 2928 | ||||

Slope of FSD | –2.15 | –2.09 | ||||

July 8 | Area (km^{2}) | 16 | 41 | 48 | 52 | 26 |

Num. of floes | 1460 | 4075 | 2322 | 2169 | 1011 | |

Slope of FSD | –2.13 | –2.13 | –2.17 | –2.15 | –2.10 |

TerraSAR-X is a German satellite carrying an X-band (9.65 GHz) synthetic aperture radar (SAR) with several different operating modes capable of imaging the Earth through clouds and darkness. We use the StripMap single-polarization (HH) images obtained from the Center for Southeastern Tropical Advanced Remote Sensing (CSTARS) at the University of Miami (https://www.cstars.miami.edu/). Images are approximately 60 × 30 km with 1.25 m pixel size. Images were acquired by CSTARS through careful tracking of four target locations marked by buoys (termed cluster #1, #2, #3, #4) as they drifted on the sea ice in the Beaufort Sea in July and August 2014 as part of the ONR MIZ Program (Lee et al., 2012). Image locations range from 73–75°N latitude and 140–160°W longitude. See Hwang et al. (2017a) for further details. We searched for TerraSAR-X images acquired within 5 days of MODIS images and within 100 km from image center to center. We found a total of 9 TerraSAR-X images that match with 3 MODIS images (July 19, July 31, and August 2, 2014). We analyzed the FSD in these 9 images; see Sections 3.3 and 4.4 below.

## 3. Methods

### 3.1. MODIS images

The steps for identifying ice floes in MODIS images are given below. Our method is similar to that of Paget et al. (2001) and Steer et al. (2008); a full discussion is given in Section 5.1.

(1) Delineate cloud-free regions manually based on visible and near-infrared bands 3-6-7, which allow discrimination of clouds over ice and snow. All subsequent processing is applied to the cloud-free regions only (see Figure 1a and b). The images used here are uncalibrated, and compositing of multiple passes can create sharp gradients. Our regions were chosen to avoid such gradients.

(2) Create a binary image to separate ice and water. Ice is relatively bright, while water is relatively dark, but choosing a threshold to apply across an entire image does not work in general because of large-scale gradients in brightness. Therefore, a local threshold must be chosen. We applied a low-pass two-dimensional filter (Gaussian-shaped kernel) with length scale 10–30 km to a cloud-free region to create a smoothed image, then subtracted the smoothed image from the original image. Positive differences were assigned a value 1 (ice); negative differences were assigned a value 0 (water). This procedure creates a binary image (see Figure 1c). The idea is that ice pixels are brighter than the average of their local neighborhood, and water pixels are darker. This method is similar to dynamic thresholding (Bernsen, 1986). We have conducted sensitivity experiments to determine how the degree of smoothing affects the final FSD; see Results below.

(3) Apply a morphological erosion operation (Serra, 1982) to the binary image in order to separate ice floes that are touching. Our standard procedure was to erode a one- or two-pixel-wide band from the binary ice regions. This procedure also has the effect of eliminating small and thin floes (e.g., 2 × *m* collections of ice pixels, for any *m*). We have conducted sensitivity experiments to see how the amount of erosion affects the final FSD; see Results below.

(4) Label the floes in the resulting eroded binary image using a recursive algorithm that identifies groups of connected pixels.

(5) Re-grow the floes outward using a variation of the morphological dilation operation (Serra, 1982) in which only those pixels that were originally assigned as ice are allowed to re-grow. This restores the floes to their original sizes (see Figure 1d).

From the re-grown, labeled image, properties of each floe are then easily calculated, such as area and caliper diameter. The caliper diameter is the distance between two parallel lines (calipers) that are just tangent to opposite sides of the floe. Averaging over all orientations of the calipers gives the mean caliper diameter (see Figure 2). We used the mean caliper diameter (MCD) as the measure of floe size, as in Rothrock and Thorndike (1984). The smallest floes we were able to identify in the MODIS images are about 2 km in size; the largest floes we were able to identify in sufficient numbers are about 30 km in size.

For any cloud-free region, the above procedure gives a set of floe sizes at one location in space and time. We applied it to 116 regions in 2013 and 140 regions in 2014. See Table 2 for statistics of the region sizes and the number of floes in each region.

### 3.2. MEDEA images

In order to identify ice floes in MEDEA images, we implemented an iterative procedure in which each iteration consists of one pass through the above steps. The first iteration identifies the largest floes, the next iteration identifies smaller floes, and so on, until the final iteration identifies the smallest floes. The iterative procedure was necessary because the very fine detail in the MEDEA images required more erosion to identify larger floes and less erosion to identify smaller floes.

The first step of the procedure, as with MODIS images, is to manually delineate a cloud-free region. The first iteration then begins by creating a binary image. A simple threshold is chosen (manually) without the smoothing step, which is not necessary because there are no large-scale gradients in brightness across the MEDEA images. The threshold is chosen by experimentation based on the visual appearance of ice floes and open water in the image.

The first erosion operation erodes a band of 40 to 100 pixels from the binary ice regions. This tends to leave a set of ice floes with many internal “holes” (water pixels) because spurious isolated water pixels within ice floes grow larger during erosion, and because there are true melt ponds (water pixels) within ice floes. Therefore the next step is to fill in the holes: water regions that are entirely surrounded by ice are changed to ice pixels. The purpose is twofold: (1) melt ponds are properly part of ice floes and therefore must be filled in; and (2) there can be ice pixels within melt ponds, and these would later be classified as ice floes themselves if they were not eliminated. Filling in the water holes prevents ice floes within ice floes. Note that filling the holes does not affect the calculation of the MCD of a floe. Filling holes was not necessary when processing the MODIS images because the amount of erosion (1 or 2 pixels) was minimal.

After the holes are filled, the floes are then labeled and re-grown using the same algorithm as for the MODIS images. The first iteration ends with the largest floes being tagged for removal in the next iteration.

The second iteration begins with the original binary image but with the tagged floes removed (they are set to water pixels). A less severe erosion operation is then used to separate touching floes, followed by hole-filling, labeling, and re-growing. The second iteration ends with the largest floes being tagged for removal in the next iteration.

The procedure typically continues for 3 or 4 iterations. The final iteration ends with all remaining floes being tagged. Thus each iteration contributes a set of tagged floes to the total.

The above method was applied to subsets of the MEDEA images (Table 3). Each subset was analyzed separately. Floe sizes identified in the subsets ranged from about 10 to 5000 m. Figure 3 illustrates the steps used in identifying floes for subset #1 of the image of July 8. After all the floes were identified, we calculated their MCDs, then analyzed the FSD.

In order to test the adequacy of the above method for producing accurate FSDs, we also manually identified 4796 floes in subset #1 of the image of July 8 (Figure 4). We then calculated the MCDs of the floes and compared the size distribution to that obtained from the method above. Results are given in Section 4.4.

### 3.3. TerraSAR-X images

Processing of TerraSAR-X images is described in detail in Hwang et al. (2017b). Subsets of size 30 × 30 km centered on reference buoys were first extracted. The resolution was reduced to 8.33 m pixel size by sub-sampling, followed by speckle filtering. They then applied the kernel graph cut algorithm to segment the image into water and ice regions, followed by the distance transform and the watershed transform to split touching ice floes, and a rule-based post-processing step for the revalidation of floe boundaries. Final results were subjected to manual checking and correction. The smallest floes that are adequately identified by the algorithm are ~100 m in diameter. We calculated the MCDs of the floes and compared the FSD from each of the 9 TerraSAR-X images with the FSD of the 3 MODIS images that are close to them in space and time. Results are given in Section 4.4.

### 3.4. Analysis of the floe size distribution

Having calculated a set of floe sizes (*x*) from an image, the next step is to characterize the distribution, *n*(*x*). Previous studies have reported power-law FSDs (Table 1), so we proceeded as follows to: (i) bin the floe sizes over an appropriate range *a* ≤ *x* ≤ *b* and plot the number of floes in each bin to check visually whether the data could plausibly be described by a power law; (ii) if a power law appears plausible, calculate the best-fitting exponent –*α* from the assumed power-law form *n*(*x*) = *cx ^{–α}* for

*a*≤

*x*≤

*b*; and (iii) apply a statistical goodness-of-fit test to determine whether the power law with best-fitting exponent does indeed describe the distribution of floe sizes sufficiently well.

Two common methods of binning are linear (all bins of equal width) and logarithmic (successive bins increase in width by a constant factor). In either case, if *x _{k}* is the center floe size of bin

*k*, and

*n*is the number of floes in bin

_{k}*k*divided by the bin width, then in a log-log plot of

*n*vs.

_{k}*x*the points would fall along a straight line if the data were power-law distributed. This visual check might then suggest using a least squares method to estimate the power-law exponent –

_{k}*α*, which would be the slope of the best-fitting line. However, that would not be the optimal method for estimating the exponent.

In the companion paper (Stern et al., 2018) we summarized results from White et al. (2008) and Clauset et al. (2009), who investigated the best way to determine the exponent of a power law from data. They found the following. (1) Using evenly spaced bins and applying an ordinary least squares method to the binned data in log-log space leads to a biased estimate of the exponent. This method should never be used. (2) Using logarithmically spaced bins and applying an ordinary least squares method to the binned data in log-log space can lead to a reasonably accurate estimate of the exponent, provided there are an adequate number of samples in each bin. (3) The best method for determining the exponent is to use the Maximum Likelihood Estimate (MLE), which does not require binning at all. According to White et al. (2008), the MLE of the power-law exponent has been shown mathematically to be the minimum variance unbiased estimator. In this study we used the MLE method for estimating power-law exponents; formulas are given in the companion paper (Stern et al., 2018, Appendices A and B).

Having calculated the best-fitting exponent –*α*, the question remains as to whether the power-law model *n*(*x*) = *cx ^{–α}* is in fact a good description of the distribution of the data. The Kolmogorov-Smirnov (KS) statistic is a measure of the distance between two probability distributions, commonly used for non-normal data (Clauset et al., 2009). It is simply the maximum distance between the cumulative distribution functions (CDFs) of the data and the model:

*D*= max|

*S*(

*x*) –

*N*(

*x*)|, where

*S*(

*x*) is the CDF of the data and

*N*(

*x*) is the CDF of the power-law model, i.e., the integral of

*n*(

*x*) with respect to

*x*. However, Clauset et al. (2009) noted that the KS statistic is relatively insensitive at the extreme limits of

*x*because the CDFs necessarily tend to zero and one there. To avoid this problem, the KS statistic can be reweighted to be uniformly sensitive across its range: $D=max|S(x)\u2212N(x)|/N(x)(1\u2212N(x))$

*S*(

*x*), and the CDF of the best-fitting power-law model,

*N*(

*x*). We denote by

*D*

_{obs}the value of

*D*calculated from the observed floe size data and the best-fitting power law model. We then used a bootstrap method (as suggested in Clauset et al., 2009) to determine whether

*D*

_{obs}is so large that we must reject the hypothesis of power-law distributed data. Let

*n*be the number of floes in the data set being tested. We run 1000 simulations, in each of which

_{f}*n*values of

_{f}*x*are drawn from a power-law distribution on [

*a, b*] with exponent –

*α*, and

*D*is calculated from this synthetic data set. The end result of the simulations is a distribution of 1000 values of

*D*, giving the sampling variability from a power-law distribution with exponent –

*α*. If the observed value

*D*

_{obs}from the actual floe size data is greater than the 95

^{th}percentile of the distribution of

*D*then we reject the hypothesis that the actual data are power-law distributed.

## 4. Results

### 4.1. Power-law fit to floe size distribution in MODIS imagery

Using the MLE method to estimate power-law exponents and the KS goodness-of-fit test described in Section 3.4, we find that only 9 out of 116 floe size data sets are rejected as not being power-law distributed in the 2013 MODIS imagery, and only 3 out of 140 floe size data sets are rejected in the 2014 MODIS imagery. We conclude that in general the FSD in the Beaufort and Chukchi seas follows a power-law distribution for floes from 2 to 30 km in size.

Figure 5 shows typical results for different dates and locations in the Beaufort Sea. The FSD in May (Figure 5a) follows a power-law distribution at both locations, with mean exponent (or slope in log-log space) –2.0. By July (Figure 5b) the FSD has changed to a power law with a steeper slope (–2.9), indicating fewer large floes relative to the number of small floes. This change occurs as the largest floes break apart due to mechanical and thermodynamic forcing from wind, currents, waves, internal stress, and melting. By October, when freeze-up begins (Figure 5c), the slope has become shallower (–2.2). This behavior is consistent with the fact that small floes have melted at a faster rate than large floes, depleting the small end of the distribution faster than the large end. In all three panels, the slopes at the two different locations are very consistent with one another.

### 4.2. Seasonal cycle of power-law exponent in MODIS imagery

How does the FSD change over the course of the year? Figures 6 and 7 show the evolution of the slope of the FSD from spring to summer to fall. It is clear that the general pattern of Figure 5 holds here: the slopes are shallow in spring and steep in summer, meaning that there are more large floes in spring and fewer large floes in summer, relative to the number of small floes, due to the break-up of the largest floes. Then, from summer to fall, the slopes become slightly shallower again, presumably because the summer melt has preferentially removed the smallest floes. This evolution is also evident in the images themselves.

Figure 8 shows the monthly mean slope of the FSD in 2013 and 2014, along with ±1 standard error of the mean. (In the calculation of the mean and standard error, each slope is weighted by the area of the cloud-free region that it represents). The seasonal cycle is evident, and it is remarkably consistent between the two years. Table 4 gives the data displayed in Figure 8. The differences between the mean May and June exponents, and between the mean June and July exponents, are statistically significant in both years, indicating that June is a transition month from spring to summer.

Year . | Quantity^{1}
. | March . | April . | May . | June . | July . | August . | Sept . | Oct . |
---|---|---|---|---|---|---|---|---|---|

2013 | Mean slope | –2.00 | –1.90 | –2.09 | –2.52 | –2.75 | –2.79 | –2.62 | –2.25 |

Std. error | 0.096 | 0.049 | 0.029 | 0.079 | 0.074 | 0.098 | 0.092 | 0.067 | |

N | 8 | 3 | 23 | 15 | 15 | 14 | 22 | 6 | |

2014 | Mean slope | –2.14 | –1.92 | –2.19 | –2.42 | –2.78 | –2.81 | –2.74 | – |

Std. error | 0.090 | 0.038 | 0.044 | 0.049 | 0.061 | 0.069 | 0.111 | – | |

N | 11 | 12 | 34 | 29 | 15 | 23 | 14 | – |

Year . | Quantity^{1}
. | March . | April . | May . | June . | July . | August . | Sept . | Oct . |
---|---|---|---|---|---|---|---|---|---|

2013 | Mean slope | –2.00 | –1.90 | –2.09 | –2.52 | –2.75 | –2.79 | –2.62 | –2.25 |

Std. error | 0.096 | 0.049 | 0.029 | 0.079 | 0.074 | 0.098 | 0.092 | 0.067 | |

N | 8 | 3 | 23 | 15 | 15 | 14 | 22 | 6 | |

2014 | Mean slope | –2.14 | –1.92 | –2.19 | –2.42 | –2.78 | –2.81 | –2.74 | – |

Std. error | 0.090 | 0.038 | 0.044 | 0.049 | 0.061 | 0.069 | 0.111 | – | |

N | 11 | 12 | 34 | 29 | 15 | 23 | 14 | – |

^{1}N is the sample size. Std. error is the standard error of the mean. Values are plotted in Figure 8.

### 4.3. Sensitivity to parameters in floe identification algorithm

We tested the sensitivity of the MODIS floe identification algorithm to two parameters: the length scale of the low-pass filter that goes into constructing the binary image (see step (2) in Methods), and the amount of erosion that is used to separate floes that may be touching (see step (3) in Methods). Changing the smoothing length scale effectively changes the size of the local neighborhood whose average brightness is subtracted from every pixel in the original image to create a difference image. Increasing the length scale creates a binary image with less small-scale structure, while decreasing the length scale creates a binary image with more small-scale structure. Changing the amount of erosion changes the ability of the algorithm to separate adjacent floes, and the ability to identify the smallest floes. The purpose of the sensitivity experiments was not to explore the entire parameter space but rather to see how small changes in the values of the parameters affected the slope of the FSD.

In the first set of experiments, we increased the length scale of the low-pass filter by 5% and found that the RMS change in the slope of the FSD was 0.8% in both 2013 and 2014. We then decreased the length scale by 5%. After removing one outlier from 2013, we found that the RMS change in the slope of the FSD was 1% in both 2013 and 2014. We conclude that the slope of the FSD is not sensitive to small changes in the length scale of the low-pass filter that goes into constructing the binary images.

In the next set of experiments, we varied the amount of erosion that precedes the labeling of distinct floes. The erosion radius specifies the size of the kernel that is used to erode ice pixels. In 48 of the 140 cloud-free regions from 2014, we used a default erosion radius of 1 pixel. For the sensitivity experiments, we changed the radius to 2 pixels. We found that the RMS change in the slope of the FSD was 8%. In the other 92 cloud-free regions from 2014, we used a default erosion radius of 2 pixels. For the sensitivity experiments, we changed the radius to 1 pixel. We found that the RMS change in the slope of the FSD was 6%. We note that in the seasonal cycle of the mean slope of the FSD (Figure 8 and Table 4), the change from April (–1.9) to August (–2.8) is 47%. We conclude that the slope of the FSD is not very sensitive to small changes in the erosion radius.

### 4.4. Extension of power law to smaller scales using high resolution imagery

The range of floe sizes detected in the MODIS images is 2–30 km, from which we have derived FSDs that follow power-law relationships. We would like to extend the FSDs to smaller floe sizes for the purpose of calibration and validation of models that include floe sizes at the meter scale (e.g., Horvart and Tziperman, 2015; Zhang et al., 2016). This extension can be achieved using satellite imagery with higher spatial resolution. Two questions arise: do the FSDs at smaller scales follow power-law relationships, and (if so) are the exponents of those power laws the same as the exponents at larger scales? If the answer to both questions is yes, then a single power law would characterize the FSD across a wide range of scales. This result would have two benefits: it would greatly simplify the implementation of the FSD in models, and it would imply that the FSD observed at scales ≥ 2 km could be extrapolated to smaller scales where it is much more difficult to observe due to sparser spatial and temporal satellite coverage.

To answer these questions we analyzed three high-resolution MEDEA images from 2014 that coincide in space and time (±3 days) with MODIS images (see Figure 9 for locations and outlines). Within each MEDEA image we identified cloud-free subsets ranging in size from 16 to 101 km^{2} (Table 3). We calculated the FSD for every subset using the procedure described in Section 3.4; we also calculated the FSD for subset #1 of July 8 based on manually identified floes. Results for the subsets of July 8 are shown in Figure 10. All the FSDs follow power-law distributions with very similar exponents, ranging from –2.10 to –2.17. This result gives us confidence that any one subset is representative of the image as a whole. The exponent of the FSD of subset #1 from manually identified floes (–2.17) is nearly the same as the exponent from the automated algorithm (–2.13), the difference being less than 2%. This result gives us confidence that the floe identification algorithm works properly.

A quantitative comparison of the power-law exponents of the FSD from near-coincident MEDEA and MODIS images (Figure 11) shows that the exponents do not agree, at least within the narrow margins of sampling variability, although in a log-log plot (Figure 12) it is difficult to visually distinguish the differences in the slopes. Note that the area covered by each MEDEA scene is very small compared to the area covered by the corresponding MODIS scene (Figure 9). The MEDEA FSD is not necessarily representative of the entire area covered by the MODIS scene; small-scale spatial variability within the MODIS scene may account for some of the difference seen in Figure 11.

We made a further comparison of MODIS-derived FSD with that from TerraSAR-X during July and August 2014. We found that the FSD from TerraSAR-X images follows a power law over floe sizes from 450 to 3000 m. (This range is similar to that of Hwang et al. (2017a), who used the same upper limit, and a lower limit varying from 250 to 980 m). Comparison of the power-law exponent of the FSD from TerraSAR-X images with near-coincident MODIS images shows excellent agreement in two cases (Figure 13), namely those with the smallest separation distance between image pairs. More inter-comparisons of low- and high-resolution imagery would be useful to clarify the conditions under which the FSD is well described by a single power law over the range of scales from 10 m to 30 km. For example, Perovich and Jones (2014) calculated that 4 m of summertime lateral melt would cause a power-law FSD to become slightly shallower for floes smaller than 100 m, giving a slight deviation from power-law behavior at the smallest scale.

## 5. Discussion

### 5.1. Identification of ice floes in images

Previous studies have used a variety of methods to identify ice floes in digital images or photographs. Rothrock and Thorndike (1984) and Lensu (1990) manually digitized floe boundaries. Banfield and Raftery (1992) created an “almost fully automatic method” involving the morphological erosion operation followed by a clustering algorithm to identify floe outlines. Soh et al. (1998), Holt and Martin (2001), and Wang et al. (2016) used a segmentation procedure followed by the morphological “restricted growing” concept, with manual correction of floe edges where they touched each other. Toyota et al. (2006, 2011, 2015) simply stated that “each sea ice floe was extracted according to its brightness,” which did not work well when ice floes touched each other, in which case the edges were manually corrected. Lu et al. (2008) used a similar method, with manual edge correction. Paget et al. (2001) and Steer et al. (2008) used morphological erosion to separate floes, followed by expansion, similar to our method in Section 3.1. Zhang and Skjetne (2015) used the “gradient vector flow snake algorithm” to identify floe boundaries based on the equilibrium of force fields constructed from the image data.

Inoue et al. (2004) adopted an entirely different approach to determining the FSD that does not require identifying ice floes at all. In a single aircraft flight of distance 360 km, they acquired 3600 digital photographs covering 80 m × 120 m each. For each photo (or frame) they applied a threshold to separate ice and water, then used an edge detection method to identify ice edge pixels and hence calculate the total floe perimeter. Assuming that a frame contains *n* circular floes of radius *r*, the total floe perimeter would be *L* = 2*n*π*r* and the ice concentration would be *IC* = *n*π*r*^{2}/*S*, where *S* is the area of the frame. Thus, by observing *L* and *IC* they were able to calculate *n* and *r* for each frame. The collection of 3600 values of *r* constituted their data on floe size, with diameters 2*r* ranging from 7 to 53 m (Table 1). From these data they constructed two FSDs, one for the internal ice zone (<250 km from coast) and one for the MIZ (>250 km from coast).

Perovich and Jones (2014) applied a similar concept that does not require identifying individual floes. First they created ten digital photo mosaics from hundreds of individual photos. In each mosaic, they calculated total floe perimeter (*L*) and ice concentration (*IC*), then assumed that the cumulative FSD followed a power law and thus solved for the power-law exponent and normalization constant in terms of the observed *L* and *IC*. In contrast to Inoue et al. (2004), this method yields a complete FSD for each mosaic, rather than just a single floe size. Each mosaic covered an area of 100–200 km^{2}, and the FSD was estimated to be valid for floe sizes ranging from 10 to 10,000 m (Table 1).

Aside from Inoue et al. (2004) and Perovich and Jones (2014), all of the methods for floe identification employ the same general strategy: separate ice from water, then identify floes either as connected sets of ice pixels or by finding their boundaries. Floe identification in images remains a complex problem. Most algorithms are not fully automated, but rather involve some interactive or manual steps, and visual checking. Our method is no exception. The new feature of our method for MEDEA images (Section 3.2) is the iterative procedure whereby large floes are identified first, then smaller and smaller floes. Fine detail is not needed to find the largest floes, but it is needed to find the smallest floes. Our iterative procedure takes advantage of that fact, eliminating the need for some manual corrections that may otherwise have been necessary.

### 5.2. Seasonal cycle of power-law exponent

Hwang et al. (2017a) and Perovich and Jones (2014) (PJ14 hereafter) are the only observational studies besides the present one in which the seasonal cycle of the power-law exponent of the FSD has been calculated and plotted. The PJ14 exponent, which is associated with the cumulative FSD, was relatively constant at about –2.0 during June and the first half of July, then dropped to about –2.2 in August before rising to –2.1 in September. Values were derived from 10 digital photo mosaics in the Beaufort Sea in 1998: 7 in June–July, 2 in August, and 1 in September. Broadly speaking, their results are consistent with the present study: the slope (in log-log space) of the FSD begins to steepen in mid-July, reaches a maximum steepness in August, and becomes slightly less steep in September. This behavior reflects floe break-up during June–July and preferential melting of small floes during August–September (Arntsen et al., 2015). Steep slopes indicate relatively more small floes than large floes.

If we attempt to make a direct quantitative comparison of PJ14 with the present study by converting our non-cumulative power-law slopes to cumulative FSDs by adding 1 (since the integral of *x ^{–α}* is proportional to

*x*

^{–α+1}), our slopes during June–September (mean values –2.4 to –2.8, Figure 8) become –1.4 to –1.8, considerably shallower than PJ14. Given the small sample size of PJ14 (only 3 values during August–September), we should not necessarily expect a close match with our mean values, which are based on sample sizes of 36 (August–September 2013) and 37 (August–September 2014; see Table 4). Other differences between PJ14 and the present study that could account for different FSDs include different methods of deriving the FSD from image data, and different years (1998 vs. 2013–2014) that span a period of substantial change in sea ice in the Beaufort Sea, which now experiences earlier spring break-up, later fall freeze-up (Stern and Laidre, 2016), and generally thinner ice than in 1998 (Lindsay and Schweiger, 2015). Wave action is also increasing during the longer, summer open-water season (Thomson et al., 2016), which affects the FSD near the ice edge, and the marginal ice zone is becoming wider in summer (Strong and Rigor, 2013; Strong et al., 2017).

As with PJ14, the time series of FSD from TerraSAR-X images (Hwang et al., 2017a) offers a Lagrangian view of the ice cover as the satellite follows clusters of buoys on the drifting sea ice. The power-law exponent of the FSD from TerraSAR-X becomes increasingly negative from July to August, 2014, indicating a shift toward more small floes relative to large floes (i.e., break-up). The mean power-law exponent for July from TerraSAR-X images is –2.57 (from 22 images), somewhat shallower than the value of –2.78 derived from MODIS images (Figure 8 and Table 4). The mean power-law exponent for August from TerraSAR-X is –3.37 (from 22 images), quite a bit steeper than the MODIS-derived value of –2.81 (Figure 8 and Table 4), possibly reflecting the fact that the TerraSAR-X images were closer to the ice edge where break-up would have been more vigorous. In general, the evolutionary trend of the power-law exponent from TerraSAR-X during July–August is consistent with this study and with PJ14.

In a study of sea-ice deformation in the Beaufort Sea, Stern and Lindsay (2009) found a power-law relationship between the mean deformation rate and the spatial scale over which it is measured, in which the power-law exponent became more negative in summer as the ice pack weakened and internal stresses were not as readily transmitted over long distances. Clearly there is a connection between sea-ice deformation and floe break-up, but the similar evolution of the power-law scaling of sea-ice deformation and the power-law scaling of the FSD over the course of spring–summer break-up has not been studied.

### 5.3. Representation of FSD and comparison with previous studies

Table 1 summarizes many of the FSD studies since the foundational work of Rothrock and Thorndike (1984). Inter-comparison of results is complicated by the fact that some researchers express the FSD as a non-cumulative density function (number density), while some use the cumulative distribution function, which is the integral of the density function. In the companion paper to this one (Stern et al., 2018), we consider these two representations of the FSD and review the results of the studies in Table 1, finding agreement between some studies and disagreement between others. The inconsistencies are possibly attributable to physical processes that affect the FSDs differently from one time and place to another, or sampling variability, or methodological differences between studies, or the inadequacy of power laws as models of the FSD. Analyzing a large number of images and reporting mean slopes (as in Figure 8 and Table 4) reduces the uncertainty due to sampling variability, and applying a goodness-of-fit test (as described in Section 3.4) quantifies the extent to which power laws are adequate models of the FSD.

We have used the non-cumulative FSD in our calculations because we find it more intuitive to understand than the cumulative FSD. Furthermore, as discussed in Stern et al. (2018), if the non-cumulative FSD follows a power law on a finite domain, then the cumulative FSD does *not* follow a power law on that domain; it has concave-down curvature in log-log space, making it difficult to uncover the underlying power-law behavior of the non-cumulative distribution from the cumulative form.

## 6. Conclusions

(1) We find that a power law describes the non-cumulative FSD in the Beaufort and Chukchi seas reasonably well, in that 244 out of 256 floe-size data sets derived from MODIS images were consistent with a power-law representation according to a Kolmogorov-Smirnov goodness-of-fit test.

(2) The monthly mean exponent (or slope in log-log space) of the FSD goes through a seasonal transformation during spring and summer (Figure 8). The slope steepens in spring as large floes break up, reducing the number of large floes relative to the number of small floes. The slope peaks in July or August and then becomes shallower, presumably as small floes preferentially disappear due to lateral melting. June is the transition month from the spring distribution to the summer distribution. The seasonal cycle is reflective of a broad area (Figures 6 and 7) and is not limited to the marginal ice zone.

(3) Our data suggest that a single power law may suffice to describe the FSD across floe sizes from 10 to 30,000 m. This is based on visual inspection of FSDs from a limited number of high-resolution images, not from a strict statistical test.

(4) Our results have been used for calibration and validation of an ice/ocean model that explicitly models the sea-ice FSD in the Arctic (Zhang et al., 2015, 2016). Two approaches to model comparison are possible: (i) observed and modeled FSDs can be compared directly (bin by bin), without regard for whether either one follows a power law or not; and (ii) observed and modeled power-law exponents can be compared. The second approach would be appropriate if the FSD in the model were parameterized as a power law *a priori*, with an evolving exponent.

## Data Accessibility Statement

This study generated data on the sea-ice floe size distributions in the Beaufort and Chukchi seas. To obtain the data reported in this paper, please contact the lead author.

## Acknowledgments

We thank NASA for the MODIS images and USGS for the MEDEA images. We thank two reviewers for comments that led to a much improved paper.

## Funding informations

This work was funded by the Office of Naval Research under the programs *Emerging Dynamics of the Marginal Ice Zone*, grant N00014-12-1-0112 (Stern, Schweiger, Stark, Zhang, Steele) and *Seasonal Ice Zone Reconnaissance Surveys*, grants N00014-12-1-0232 and N00014-17-1-2545 (Schweiger, Stark, Steele), and by the NASA Cryospheric Sciences Program, grant NNX17AD27G (Stern, Schweiger, Zhang).

## Competing interests

The authors have no competing interests to declare.

## Author contributions

Contributed to conception and design: All

Contributed to acquisition of data: Schweiger, Stark, Hwang

Contributed to analysis and interpretation of data: All

Drafted and/or revised the article: Stern, Schweiger, Zhang, Steele

Approved the submitted version for publication: All