the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
Estimating the sea ice floe size distribution using satellite altimetry: theory, climatology, and model comparison
Lettie A. Roach
Rachel Tilling
Cecilia M. Bitz
Baylor FoxKemper
Colin Guider
Kaitlin Hill
Andy Ridout
Andrew Shepherd
In seaicecovered areas, the sea ice floe size distribution (FSD) plays an important role in many processes affecting the coupled sea–ice–ocean–atmosphere system. Observations of the FSD are sparse – traditionally taken via a painstaking analysis of ice surface photography – and the seasonal and interannual evolution of floe size regionally and globally is largely unknown. Frequently, measured FSDs are assessed using a single number, the scaling exponent of the closest powerlaw fit to the observed floe size data, although in the absence of adequate datasets there have been limited tests of this “powerlaw hypothesis”. Here we derive and explain a mathematical technique for deriving statistics of the sea ice FSD from polarorbiting altimeters, satellites with subdaily return times to polar regions with high alongtrack resolutions. Applied to the CryoSat2 radar altimetric record, covering the period from 2010 to 2018, and incorporating 11 million individual floe samples, we produce the first panArctic climatology and seasonal cycle of sea ice floe size statistics. We then perform the first panArctic test of the powerlaw hypothesis, finding limited support in the range of floe sizes typically analyzed in photographic observational studies. We compare the seasonal variability in observed floe size to fully coupled climate model simulations including a prognostic floe size and thickness distribution and coupled wave model, finding good agreement in regions where modeled ocean surface waves cause sea ice fracture.
 Article
(2900 KB)  Fulltext XML

Supplement
(108 KB)  BibTeX
 EndNote
Earth's polar oceans are covered with sea ice: a thin, heterogeneous interface that plays an important role in the coupling between ocean and atmosphere. Sea ice is a collection of many individual pieces, called floes, which may be characterized in terms of a horizontal length scale, their “size”. On the large scales relevant to global climate modeling, the statistical variability of floe size is described using the floe size distribution (FSD; Rothrock and Thorndike, 1984).
The FSD is an important property of the sea ice cover that influences the multiscale temporal and geographic variability of sea ice, akin to the grain size in sedimentology or particle size distribution in atmospheric chemistry. The scale of individual floes plays a role in many seaicerelated processes: sea ice melt rate (Steele, 1992; Horvat et al., 2016; Horvat and Tziperman, 2018), the evolution of the oceanic mixed layer (Manucharyan and Thompson, 2017), atmospheric boundary layer exchange (Birnbaum and Lüpkes, 2002; Lüpkes and Birnbaum, 2005; Tsamados et al., 2014), the sea ice response to applied stress (Feltham, 2008; Wilchinsky and Feltham, 2011), and the propagation of waves into the ice (Squire et al., 1995; Squire, 2007; Smith and Thomson, 2016). The importance of the sea ice FSD has led to the development of diagnostic FSD models of varying complexity (Williams et al., 2013; Zhang et al., 2015; Bateson et al., 2019) and a prognostic floe size and thickness distribution (FSTD) scheme (Horvat and Tziperman, 2015; Roach et al., 2018a).
Despite the potential relevance of sea ice floe size to polar climate evolution, there remain no climatescale assessments of average floe size or the FSD. The observational record of floe statistics derives from visual imagery localized in space and time (i.e., Rothrock and Thorndike, 1984; Toyota et al., 2006, 2011; Steer et al., 2008) or from repeat measurements in the same region over multiple months (Hwang et al., 2017; Stern et al., 2018a), which may subsequently be used to compile a seasonal cycle of the FSD (Perovich and Jones, 2014; Stern et al., 2018a). FSD measurements are obtained by identifying individual floes within a twodimensional image of the sea ice surface. Because floe sizes span several orders of magnitude, accurate representations of the FSD – even in relatively small geographical domains and in perfect lighting and surface conditions – require high resolution and high observational coverage. Nearly all measurements of the FSD have been made in accordance with a “powerlaw” scaling hypothesis commonly used to describe multiscale systems (Mandelbrot and Wheeler, 1983), in which the resulting FSD is fit to a straight line in logarithmic coordinates, whose slope, α, is reported as an intrinsic property of the floe mosaic. There is large uncertainty in these scaling coefficients, the range they apply over, and their applicability and origin (Herman, 2011; Horvat and Tziperman, 2017; Herman et al., 2018; Stern et al., 2018b). Improvements in the quality and quantity of available FSD data are needed before arriving at consensusderived FSD statistics to guide and assess model performance.
Here we outline a method that exploits satellite radar altimetry to construct the FSD and its moments across polar regions with subkilometer spatial resolution, subdaily temporal resolution, and spanning multiple orders of magnitude in size. Altimeters, like the ones carried on the Envisat, ICESat, CryoSat2, and ICESat2 satellites, make repeated, frequent passes over polar oceans, and substantial efforts have been made to process the satellite returns to discriminate between open water, floes, and leads. The altimetric returns have found many uses, including reconstructing the sea ice thickness field (Laxon et al., 2013; Tilling et al., 2016, 2018) and ocean surface circulation under sea ice (Peacock and Laxon, 2004; Armitage et al., 2018). Fields inferred from altimetry have led to advances in understanding polar systems: from forecast and climate prediction (Day et al., 2014) to model validation (Schröder et al., 2018; Allard et al., 2018) to climate change studies (Laxon et al., 2003; Kwok, 2018), and have been evaluated and validated using field campaign data (Skourup et al., 2017; Sandberg Sorensen et al., 2018; Tilling et al., 2018).
Onedimensional measurements of sea ice properties, like alongtrack altimetric measurements of ice open water, have long been sought to describe the twodimensional ice surface. Rothrock and Thorndike (1984) originally described a method for reconstructing the sea ice floe size distribution in a region using straightline measurements over the geometry of floes. Lindsay and Rothrock (1995) later compiled the statistics of lead and ice spacings in twodimensional imagery. Other work has taken place to derive and understand the width distribution of individual leads in visual imagery and altimetry (Wadhams et al., 1988; Key and Peckham, 1991; Key, 1993; Wernecke and Kaleschke, 2015), which can be used to estimate heat fluxes and turbulent transfer between the ocean and atmosphere. To date, however, these studies have not been designed to facilitate a comparison with model data, nor have altimetric studies been used to compile floe size statistics. These objectives are the focus of this work.
We outline the mathematical theory that allows for comparison of altimetric datasets and the FSD in Sect. 2. In Sect. 3 we apply this method to a new dataset of segmented CryoSat2 seaicetype data from 2010 to 2018. Using these data we produce the first climatological maps of mean sea ice floe size and fragmentation for the Arctic Ocean. We then test the powerlaw hypothesis, finding limited support for powerlaw scaling across most of the dataset in Sect. 4. One of the key aims of the paper is to develop floe size distribution measurements that are useful for model validation and calibration. In Sect. 5, we show a proof of concept, demonstrating how altimetric data can be used to constrain and evaluate new models of the FSD, comparing the CryoSat2 FSD data to a climate model simulation with a prognostic FSTD model. We conclude in Sect. 6.
For an individual pass over sea ice by a polarorbiting satellite altimeter, return waveforms along the satellite orbit track are assigned a surface type depending on the waveform shape and coincident sea ice concentration (Tilling et al., 2018). A “floe chord” of length D is a continuous series of points identified as sea ice, covering a geographic distance D (Tilling et al., 2019a, b). Define a floe's size, r, as its “effective radius” – the square root of the floe's area divided by π (Rothrock and Thorndike, 1984; Horvat and Tziperman, 2015) We use radius instead of diameter, as appears in some other observational studies, for comparison with model output in Sect. 5. Because the satellite path is at an unknown angle with respect to the (also unknown) floe geometry, any individual floe chord measurement is not a floe size measurement. Converting between suitably processed altimetric floe chord measurements and floe size statistics is therefore the subject of this section. Details on the processing of the CryoSat2 waveform, used to produce a dataset of floe chords spanning the period 2010–2018, are outlined in Sect. 3 and Tilling et al. (2019b).
For a domain of horizontal area A, and over a period of time ΔT that corresponds to several repeat satellite passes, we bin the set of recorded floe chords to form a probability distribution S(D), which we term the “floe chord distribution” (FCD), where S(D)dD is equal to the number fraction of floe chords in A over ΔT with length between D and D+dD, and is normalized to 1. To collapse all measured chords onto a single independent scalar coordinate (D), we follow the example of turbulence statistics (Batchelor, 1953) and assume that the floe chord distribution data are homogeneous, isotropic, and stationary within the region and time data are collected. In the same region, we define the (noncumulative) number FSD P(r), where P(r)dr is the fractional number of floes with a size between r and r+dr in A, and is also normalized to 1. The FSD inherits the assumptions of homogeneity, isotropy, and stationarity from the FCD. Our objective is to relate the FCD, S(D), or quantities derived from the FCD, to the statistics of the FSD, P(r).
Bayes' theorem relates S(D) and P(r) through conditional probabilities,
The conditional probability F(r;D) relates given chord lengths to the floe size distribution that could generate them: F(r;D)dr is the probability that floes with size in the range from r to r+dR were sampled given a chord of length D. The conditional probability $\stackrel{\mathrm{\u0303}}{F}(D;r)$ relates given floe sizes to the chord length distribution they generate: $\stackrel{\mathrm{\u0303}}{F}(D;r)\mathrm{d}D$ is the probability of measuring a floe chord of length from D to D+dD given that a floe of size r was measured.
This second probability distribution $\stackrel{\mathrm{\u0303}}{F}(D;r)$ can be derived from first principles under a single assumption: that the chord length distribution that would be sampled from a set of floes of size r is independent of r (equivalently, the floe shape distribution is scaleinvariant). Formally, this requirement is
where $G\left(\mathit{\xi}\right)=G\left(\frac{D}{\mathrm{2}r}\right)$ is an unknown function that integrates to 1 over the interval from ξ=0 to 1. Under this assumption, the distribution of possible chord lengths measured from floes of size r has the same functional form independent of r. The probability distribution F(D;r) may be derived by considering the geometric relationship between straightline satellite passes and the geometry of the floes they pass over. Individual floe shapes are highly variable: making an assumption about the distribution of floe shapes may introduce biases in the statistics derived from the FCD. Yet as we prove in Appendix A, the ability to derive FSD statistics from the FCD does not depend on the precise form of $\stackrel{\mathrm{\u0303}}{F}(D;r)$ so long as the homogeneous, isotropic, stationary, and scaleinvariance assumptions are retained, and the evaluation of powerlaw scaling is in fact independent of $\stackrel{\mathrm{\u0303}}{F}(D;r)$.
To proceed and arrive at a concrete (although not general) realization of these functions, we will assume all floes are perfect circles. In assessments of the relationship between major and minor axes of individual floes, the “roundness” parameter for a floe is typically within 15 % of 1 (Rothrock and Thorndike, 1984; Toyota et al., 2011; Perovich and Jones, 2014; Gherardi and Lagomarsino, 2015; Alberello et al., 2019), suggesting that this circular assumption, while simplistic, is broadly appropriate. Nevertheless, it will likely be necessary to amend the analysis below in the future to account for more realistic shape distributions and geometries (e.g., diamonds; Wilchinsky and Feltham, 2006), regional differences in floe shape properties (such as in regions where shear stress determines fracture patterns and floe shapes; Schulson and Hibler, 1991), or to evaluate the sensitivity of the results that follow to the assumed shape distribution. Solving for $\stackrel{\mathrm{\u0303}}{F}(D;r)$ is a geometric problem that relates the possible measured chord lengths to the underlying floe size, and we solve this explicitly for circular floes here. Similar geometric problems have been identified and solved in other fields (e.g., Pons et al., 2006; Nere et al., 2007), and we therefore leave refinement of $\stackrel{\mathrm{\u0303}}{F}(D;r)$ to future work.
Consider the special case that all floes are perfect circles, illustrated in Fig. 1. Because there is no correlation between the statistics of local sea ice deformation and predetermined satellite tracks, an individual recorded floe chord, D, originating from a floe of radius r, was obtained from a satellite trajectory that crosses the floe at a random interior angle θ; thus the distribution of θ is uniform. Because of rotational symmetry, we need only consider $\mathit{\theta}\in [\mathrm{0},\mathit{\pi})$, sampled according to a probability distribution $T(\mathit{\theta};r)={\mathit{\pi}}^{\mathrm{1}}$. The length D is thus a chord of this circular floe, with $D=\mathrm{2}r\mathrm{sin}\left(\mathit{\theta}/\mathrm{2}\right)$. Accordingly,
which is a probability function that meets the above criterion (Eq. 2).
The nth moment of the floe chord distribution S(D) is defined
For any function $\stackrel{\mathrm{\u0303}}{F}(D;r)$ satisfying the scale invariance above, the righthand side may be expressed in terms of moments of P(r) (see Appendix A). For circular floes, using Eq. (3),
where $\frac{D}{\mathrm{2}r}\equiv \mathit{\xi}=\mathrm{sin}\left(x\right)$, 〈r^{n}〉 is the nth moment of P(r), and the coefficient A_{n} is
where B is the beta function. For n=0, 1, 2, or 3, then A_{n} is 1, $\frac{\mathrm{4}}{\mathit{\pi}}$, 2, or $\frac{\mathrm{32}}{\mathrm{3}\mathit{\pi}}$. Two important FSDderived quantities are derived from ratios of FSD moments, and therefore can be obtained from the FCD directly: the “representative radius” (Horvat and Tziperman, 2017; Roach et al., 2018a),
and the floe perimeter per ice area, a measure of sea ice fragmentation,
These derived quantities are useful because they require no further information about the sea ice (such as its concentration) to compare against modeled FSDs. However, both $\stackrel{\mathrm{\u203e}}{r}$ and 𝒫 can represent only those floes whose size is larger than ${r}_{\mathrm{min}}={D}_{\mathrm{min}}/\mathrm{2}$, the smallest possible floe size sampled. For perfect powerlaw distributions beginning at a scale of r_{min} or before, both metrics are functions of r_{min}. However, for the real FCDs measured here, a maximum floe size exists, and a powerlaw scaling is not found approaching r_{min}, so the use of such metrics is justified (see Sect. 4). Because of the finite sampling resolution of the altimeter, chords that would originate from floes with a diameter near the sampling resolution may not be observed, and thus 〈D^{n}〉≤A_{n}〈r^{n}〉. We explore this uncertainty in Appendix B. For a known floe size distribution, the error decreases exponentially as a function of the distributional moment being considered, though it can be large (20 % or more) in pathological cases. For distributional tails characterized by observed scaling exponents (Stern et al., 2018b), and for moments considered here, this uncertainty can be determined systematically and vanishes for measurement spacings smaller than the radius of the most common floe size. This resolution error does not affect the analysis of the powerlaw hypothesis, as that analysis is focused on the distributional tail. However, because 𝒫 is proportional to a negative moment of the FCD, it is sensitive to changes in the number of small chord lengths. Because of the measurement uncertainty for smaller chord lengths we will focus instead on $\stackrel{\mathrm{\u203e}}{r}$, which is a positive moment of the FCD.
2.1 Evaluating the floe size powerlaw hypothesis with floe chord data
Suppose the FSD P(r) has a powerlaw tail that begins at some specified value r_{1}. Then for r>r_{1}, $P\left(r\right)\equiv P(r;\mathit{\alpha},C)=C{r}^{\mathit{\alpha}}$, for an unknown coefficient C and powerlaw slope α. Integrating Eq. (1) over all r,
where the integral of the lefthand side of Eq. (1) is equal to S(D) as $\int F(r;D)\mathrm{d}r=\mathrm{1}$. Under the assumption of Eq. (2), if P is a power law, so is S(D) (Appendix A). For circular floes,
Because of the sampling resolution of the altimeter there is a minimum resolved chord scale D_{min}. If ${D}_{\mathrm{min}}\ll {D}^{*}\equiv \mathrm{2}\cdot {r}_{\mathrm{1}}$, there is an explicit solution for S(D), a powerlaw distribution over the range $({D}^{*},\mathrm{\infty})$,
where B is the beta function. The coefficient C_{α} is a multiplicative factor independent of size, and the powerlaw exponent for a FCD is the same as the exponent for FSD, where the two are related by Eq. (1).
Moments of a powerlaw tail can be evaluated explicitly (for $\mathit{\alpha}>n+\mathrm{1}$),
Then for both the FCD and FSD, the ratio of two moments is independent of the unknown coefficient C, i.e.,
valid for $n+\mathit{\u03f5}<\mathit{\alpha}$. The powerlaw coefficient can be obtained for any n, ϵ as
In the analysis below we will arbitrarily select only $n=\mathrm{0.5},\mathit{\u03f5}=\mathrm{1}$ for comparison (for scaling coefficients α>1.5, the bulk of reported powerlaw coefficients are in this range; Stern et al., 2018b). Because the observations will not be perfect powerlaw distributions, we will use ${\mathit{\alpha}}_{\mathrm{0.5},\mathrm{1}}\equiv {\mathit{\alpha}}^{*}$ as an estimator. A second estimate of the powerlaw scaling coefficient, $\widehat{\mathit{\alpha}}$, is computed via the maximum likelihood estimator (Muniruzzaman, 1957; Clauset et al., 2009; Virkar and Clauset, 2014) (details in Appendix C) as
where N is the number of chords. If the powerlaw hypothesis holds, then the two estimates of α agree, although the agreement of $\widehat{\mathit{\alpha}}$ and α_{n,ϵ} is not sufficient to confirm the powerlaw hypothesis. In the Supplement Sect. S1, we compare these two estimates when they are evaluated against synthetic datasets drawn from a true powerlaw distribution. The two agree even when the size of the data is relatively small (N<25). While in practice Eq. (13) is easy to apply, it only holds when ${\mathit{\alpha}}_{n,\mathit{\u03f5}}>n+\mathrm{1}$, and unlike the method of Clauset et al. (2009), it does not allow for a robust statistical analysis of the powerlaw fit, and should only be used when the data are assumed to follow a powerlaw already.
We apply the analytic technique described in Sect. 2 to a floe chord dataset constructed from the CryoSat2 radar altimeter processed by the Center for Polar Observation and Modelling (CPOM) over the period from October 2010 to present (CPOM data products are available at http://www.cpom.ucl.ac.uk/csopr/seaice.html, last access: 9 January 2019). CryoSat2 radar echo returns are defined as “lead”, “floe”, “open ocean”, or “ambiguous” according to waveform shape and sea ice concentration (Tilling et al., 2016, 2018), at an approximately constant alongtrack spacing D_{min}=300 m. Floe chords are defined as a continuous sequence of one or more “floe echoes”, with a gap of one ambiguous echo permitted within a floe sequence to allow for anomalous returns. A chord length is taken from the midpoint of the first to the midpoint of the last radar echo. Individual chord lengths can be underestimated when continuous floes are separated artificially by producing two or more ambiguous echoes in sequence, or when highly reflective leads dominate the waveform return close to the floe edge and cause measurement dropout (Tilling et al., 2019b). Lead contamination, or “snagging” (Armitage and Davidson, 2014), is more likely when the altimeter cuts off a small section of a floe, i.e., for small values of θ. Overestimates of chord length can also occur when ice floes are in close contact with neighboring floes. Therefore, floe chord lengths should be considered a satellitederived product, not a true measurement of floe size. The minimum chord length retrieval D_{min} is limited to the CryoSat2 footprint (∼300 m alongtrack) (see the discussion in Appendix B). However, surface discrimination via altimetry is highly accurate in months without melt ponds (Peacock and Laxon, 2004; Guerreiro et al., 2017; Quartly et al., 2019), giving confidence that floe echoes represent a coherent length of ice. More details on the details of chord identification may be found in Tilling et al. (2019b). Indeed, these raw floe chord data have been used successfully to reduce biases in altimeterobserved satellite sea ice thickness estimates from altimeters with different footprint sizes (Tilling et al., 2019b). Here we analyze the sea ice floe size distribution using that floe chord product.
Figure 2 shows an example of floe chord data for a single CryoSat2 track over the Arctic on 21 January 2014. Freeboard values for echoes discriminated as floe are plotted in Fig. 2b as a function of the alongtrack distance in kilometers, and correspond to the blue circle in Fig. 2a. Floe chords are identified as black segments in Fig. 2b. The histogram of all 741 identified chords for this single satellite pass is shown in loglog space in Fig. 2d.
The full CryoSat2 dataset examined here spans the time period from October 2010 to November 2018, and floe chords measured using the above technique are binned into the CICE sea ice model's twodimensional sea ice grid for each month and year to facilitate comparison with model products. This implies that we invoke the principles of isotropy, homogeneity, and stationarity of the FCD, required to produce such a distribution, on the length scale of the CICE model grid (O(25 km)) and timescale of a month. For every grid cell i, month m, and year y, we have a vector of floe chords $\mathit{\left\{}{D}_{i,m,y}\mathit{\right\}}$ from which we build a FCD. The base 10 logarithm of the total number of floe chords recorded in each grid cell per month is shown in Fig. 2a. Because the satellite passes are densest near the pole, the measurement density is highest near the pole as well. Figure 2c shows the number of Arctic measurements in each month. Sea ice type from CryoSat2 is not available during summer months, as melt ponds make it difficult to discriminate between leads and ponded floe surfaces, and we do not include measurements from May to September. Across the entire set of satellite tracks included here, 11 million chord lengths are recorded in the Arctic.
Figure 3a shows the seasonal cycle of Arctic representative radius over the CryoSat2 period obtained by applying Eq. (6) to the binned CryoSat2 floe chord product. Individual years are plotted as thin lines, and the climatological average is shown in red. Details on how temporal and spatial average statistics are computed are included in Appendix D. During the months of October–December, the climatological representative radius is roughly 35 % larger (7.06 km vs. 5.18 km) than during February–April. This seasonal cycle is broadly consistent across years. A possible interpretation of this seasonal cycle is that large firstyear ice pans form in October and are later fractured into smaller floes throughout the winter months. This concept is supported by observations that largescale fracturing of sea ice in the Beaufort Sea is dominated by coastal processes and therefore can only occur once sea ice freezes to the coast in midwinter (RichterMenge, 2002), although such an interpretation is speculative and must be evaluated further as this method is refined. Figure 3b shows annualaverage representative radius in red for each full year from 2011 to 2017, with thin lines corresponding to the individual months within that year. Seasonal variability is significantly larger than interannual variability. There is no statistically significant linear trend at the p=0.05 level.
The geographic variability of representative radius over the “early winter” (October–December) and “late winter” (February–April) periods is shown in Fig. 3c–d, for all grid areas. We display only those areas with at least 25 recorded floe lengths in each month during the averaging period. In Sect. S2 and Fig. S1, we examine the sensitivity of bulk FSD statistics to this threshold, finding similar seasonal cycles and climatologies. The largest representative radii in the Arctic lie in the interior Arctic near the pole, with a tongue of large floes that extends along the Canadian Arctic in late winter. There is a notable increase in representative radius with latitude. In Fig. S2, we show that this relationship cannot be explained as a result of the increasing density of measurements near the pole and may therefore be a geophysical signal. The smallest representative radii (below 1 km) lie in the Bering Strait and the Russian Arctic in early winter and in the Laptev Sea in late winter. The difference in representative radius between fall and spring is accounted for by the reduction of floe sizes in regions near the Arctic interior (see Fig. 6).
Given a collection of chord lengths, we would like to examine whether it is distributed according to a power law. Under the assumptions of Sect. 2, the scaling behavior of the FSD is the same as of the FCD (see Appendix A). We use the statistical methodology outlined in Clauset et al. (2007, 2009) and Virkar and Clauset (2014) (which we term the MLE method) to evaluate shape parameters of the most likely powerlaw fit and to test its plausibility. This method has been used to evaluate powerlaw behavior in a recent FSD model (Horvat and Tziperman, 2017) and observational studies (Hwang et al., 2017; Stern et al., 2018b) and proceeds as follows.

Lowertruncate the FCD. First identify a minimum chord scale, D^{*}, above which we hypothesize a powerlaw tail, and analyze only those floe chord measurements. We either (a) choose D^{*} to be 900 m (to reduce the impact of smallsize sampling errors discussed in Sect. 2) or (b) use the scheme described in Clauset et al. (2007) to evaluate the most likely value of D^{*} for a powerlaw tail. The length of this lowertruncated distribution is N. In the descriptions that follow, we use the subscript “all” to describe case (a) and “tail” to describe case (b).

Compute powerlaw scaling estimates and parameter uncertainty. We obtain two estimates of the FCD scaling estimate, either computing α^{*} via Eq. (13) or computing $\widehat{\mathit{\alpha}}$, and uncertainty estimates in both $\widehat{\mathit{\alpha}}$ and D^{*} via the MLE method (Eq. 14). That the two estimates of α agree is a necessary condition for the FCD (and thus FSD) to be powerlaw distributed.

Examine the plausibility of the powerlaw fit. We generate M FCDs of size N (the same number of synthetic chords as observed chords), with each synthetic FCD drawn from the hypothesized powerlaw distribution $P(\widehat{\mathit{\alpha}},{D}^{*})$. For each of these synthetic FCDs, we compute the Kolmogorov–Smirnov distance between it and the hypothesized powerlaw model that generated it, $P(\widehat{\mathit{\alpha}},{D}^{*})$. We also compute the distance between the observed FCD and $P(\widehat{\mathit{\alpha}},{D}^{*})$. A p value, p, is equal to the fraction of those M synthetic FCDs that are “further away” from the hypothesized powerlaw model than is the observed FCD. We use M=10 000, which permits computation of p within 0.005 (Clauset et al., 2009), and rule out the powerlaw hypothesis under the condition p<0.1 (Virkar and Clauset, 2014).
We note that a power law describes the scaling of a distribution's tail. Previous observational studies have discussed “double power laws” (i.e., Toyota et al., 2011), i.e., two powerlaw distributions of a different exponent joined at a specified scale. The methods employed here would capably capture the largesize powerlaw scaling but not the smallscale scaling. Such double power laws are necessarily scalevariant, and require at least three parameters to describe. The conceptual and mathematical simplicity of the powerlaw hypothesis does not apply in such a case, and we do not consider them here.
The MLE method is a rigorous test of the powerlaw hypothesis that eliminates potential human bias when interpreting observational data. To illustrate why this is important, we first consider the entire set of 11 million chord lengths recorded in the Arctic in all months (October–April), spanning a length range from 300 m to 100 km. The histogram of these floe chords is the black line in Fig. 4a (hashes on black line are the logarithmically spaced bin centers). Beginning from ${D}^{*}={D}_{\mathrm{min}}=\mathrm{900}$ m, ${\widehat{\mathit{\alpha}}}_{\mathrm{all}}=\mathrm{1.97}$ (blue line) and ${\mathit{\alpha}}_{\mathrm{all}}^{*}=\mathrm{2.05}$ (not shown). The observations are further away from synthetic data drawn from $P({\widehat{\mathit{\alpha}}}_{\mathrm{all}},{D}^{*})$ in each of the M=1000 random draws (${p}_{\mathrm{all}}=\mathrm{0}/\mathrm{1000}$) and we reject the powerlaw hypothesis for these measurements. We note that if the resolution bias explored in Appendix B proves to be larger than expected, the underrepresentation of small floe lengths may affect the analysis of the full distribution.
Examining the tail of the distribution in Fig. 4a, the maximum likelihood estimate of D^{*} is ≈15.0 km (red vertical line, vertical shaded region is the range of uncertainty for D^{*}), above which there are ∼40 000 chord length measurements between 24.7 and 99 km (0.4 % of the dataset). On the truncated FCD, ${\widehat{\mathit{\alpha}}}_{\mathrm{tail}}=\mathrm{4.65}$ (red line, dashed lines are uncertainty ranges for ${\widehat{\mathit{\alpha}}}_{\mathrm{tail}}$) and ${\mathit{\alpha}}_{\mathrm{tail}}^{*}=\mathrm{4.67}$ (not shown), similar to the largescale rolloff reported in observations (Toyota et al., 2016). Even when restricted to the FCD tail, ${p}_{\mathrm{tail}}=\mathrm{0}/\mathrm{1000}$.
Finding no statistical basis for a powerlaw fit to the tail in Fig. 4a underscores the challenge in using the human eye to observe powerlaw scaling. While the black and red lines in Fig. 4a appear similar across much of the range of sizes above 24.7 km, examining the misfit between the powerlaw estimates and the data shows that the two curves in fact differ significantly across the entire fit range. A misfit error can be defined as
where x_{i} represents the bin locations, angle brackets denote an average over the relevant bins, and P(x_{i}) represents the observed histogram values. Over the range from 24.7 to 100 km, the misfit error is 33 %. The visual agreement, misfit error, and apparent slope and shape of the distribution depend sensitively on the bin spacing and the logarithmic plotting.
Sea ice parameterizations that assume a powerlaw distribution may significantly bias sea ice statistics. The imposition of any fixed distributional shape, when FSD dynamics are scalevariant, leads to implicit nonlocal redistribution of sea ice between floe size categories (Horvat and Tziperman, 2017). To see this in practice we compare the difference in Arcticwide representative radius, $\stackrel{\mathrm{\u203e}}{r}$, which is used in parameterizations of wave attenuation and ice thermodynamics, between the mostlikely powerlaw fit to the data and the “true” value obtained via Eq. (6). The observations yield $\stackrel{\mathrm{\u203e}}{r}=$10.2 km, versus 34.5 km for the powerlaw fit. Examining only the tail of the distribution (chord lengths above 24.7 km) yields better agreement: 23.7 km for the observations and 24.4 km for the fit line. Yet this tail constitutes just 1 % of all measured chord lengths, corresponding to just 18 % of total ice area and 4.5 % of the perimeter per square meter (Eq. 7).
Segmenting the chord length data into individual months in the Arctic, there are none where p_{all}>0. Examining only the tail of each month's distribution, p_{tail}<0.1 in all months. Only in April is there a nonzero p_{tail}=0.04, for which the analysis of Fig. 4a is repeated as Fig. 4b. In April, ${\widehat{\mathit{\alpha}}}_{\mathrm{all}}=\mathrm{1.99}$, ${\widehat{\mathit{\alpha}}}_{\mathrm{tail}}=\mathrm{5.70}$, and ${D}^{*}=\mathrm{30.7}$ km. The tail consists of 1618 measured chord lengths up to 97.5 km, accounting for 8 % of the total floe area and 1.4 % of the perimeter per square meter. The misfit error between the April FCD tail and $P({\widehat{\mathit{\alpha}}}_{\mathrm{tail}},{D}^{*})$ is 76 %. Accumulating all measured chord lengths from October to May into the CESM model grid, we find zero locations that support a powerlaw distribution across the range of measurements (i.e., p_{all}>0.1). For grid areas with N>1000, we show the value of D^{*} computed using the local FCD in Fig. 4c. Values of D^{*} range from 2 km along the Russian Arctic to more than 10 km near the North Pole.
While most of the Arctic has at least 1000 total measurements across all years, FCD tails ($D>{D}^{*}$) are not as wellsampled. We investigate these tails including regions with at least 200 measured floe chords larger than D^{*}. The percentage of geographic areas with at least 1000 total measurements that have a tail with at least 200 measurements is 44 %; on average D^{*} is 5.4 km for these regions. For most of these regions we can not rule out a powerlaw tail. For the subset of regions with 1000 total measurements, 200 measurements in the tail, and where the powerlaw hypothesis cannot be ruled out, the average D^{*} is 6.5 km and average ${\widehat{\mathit{\alpha}}}_{\mathrm{tail}}$ is 3.34, within the typical range of Arctic FSD measurements (Stern et al., 2018b). In Fig. 4d we show the values of ${\widehat{\mathit{\alpha}}}_{\mathrm{tail}}$ at these locations. Colored cells are those with p>0.1 and a tail with at least 200 measurements. In Fig. 4e we show the base 10 logarithm of the MLE tail for all geographic locations. Those regions for which a power law cannot be ruled out are generally those with the largest floes and the highest sampling, clustered near the central Arctic. The weakest support for a powerlaw tail is in the Chukchi and Beaufort seas, where powerlaw floe size distributions have often been reported. We note that our choice of tail length plays an important role in whether the powerlaw hypothesis is rejected in the tail across the Arctic. For example, the fraction of Arctic regions with at least 1000 total measurements, a tail of at least 100, 200, and 400 measurements, and that does not reject the powerlaw hypothesis is 72 %, 52 %, and 15 %, respectively. The bettersampled the FCD–FSD, the more likely the powerlaw hypothesis is rejected.
Scaling coefficients can provide useful information about the distributional shape. In Fig. 5a–d we show the seasonal and interannual variability of powerlaw estimates in the Arctic. Figure 5a plots the climatology of the powerlaw scaling estimates when including all measured chord lengths in dark red (using Eq. 13) or blue (using 14). Individual years are thin red or blue lines. The two estimates disagree. Because agreement between the two estimates is necessary for the powerlaw hypothesis to be true (see Sect. 4, Sect. S1), this alone is sufficient to rule it out. There is a seasonal cycle in the powerlaw fitting to the full distribution, with α_{all} increasing (steepening) from September to January and remaining flat until April and no significant linear trend at the p=0.05 level for the annualaverage value of α_{all}. Figure 5c–d repeat this analysis on the tail of those monthly distributions. In this case, the two estimates agree well. There is a different seasonal cycle in the steepness of the distributional tail: shallowest in early winter and steeper in late winter. This indicates that the changes across the winter months may be due to a reduction of the largest floes and a steepening of the distributional tail, although there is significant interannual variability among these estimates. A similar seasonal cycle to that found in Fig. 6a, c, with an FSD that steepens from September to April, was found in image analysis of floes in the Beaufort and Chukchi seas (Stern et al., 2018b), with α≈2.5, although the distribution steepened monotonically over that period. There is no significant linear trend at the p=0.05 level in the annualaveraged FSD tail slope (Fig. 5d).
A principal aim of this work is to allow model–data comparisons and facilitate testing rapidly developing FSD–FSTD models. Here we demonstrate how such a comparison can be made and provide useful information to modelers, even in the presence of the high uncertainties in this nascent FSD reconstruction technique. With the gridded data provided above, we may now directly compare developmentstage sea ice models that incorporate FSD effects to observations. To do so, we use the Roach et al. (2018a) prognostic model for the FSD–FSTD, based on the Horvat and Tziperman (2015) theoretical FSTD framework, implemented into the CICE 5.1.2 (Hunke et al., 2015) sea ice model. The FSTD is a sea ice state variable, subject to interaction of five key physical processes: lateral growth, lateral melt, fracture by ocean surface waves, welding of floes in freezing conditions, and wavedependent new ice growth (Horvat and Tziperman, 2015, 2017; Roach et al., 2018a, b). Previously published model runs (Roach et al., 2018a) focused on the impact of the FSD on lateral melt, which is largely driven by small floes (Steele, 1992), and so floe sizes above 1 km were not considered. As a larger range of scales is resolved in the CryoSat2 observational product, we conducted a model run that extended the floe size categories to scales larger than 1 km, using 24 logarithmically spaced floe size categories from 0.5 m to 33 km.
This FSTD model simulation is coupled to a slab ocean model and the WAVEWATCH III ocean surface wave model (Tolman, 2009), forced by the JRA55 atmospheric reanalysis (JRA55, 2013) over the period from 2000 to 2016. These wavecoupled runs are branched at year 2000 from a standalone sea ice run from 1975 to 2000, spun up using repeated 1975 atmospheric forcing. Additional model physics beyond those processes outlined in Roach et al. (2018a) have been added to determine the initial size of newly formed sea ice floes as a function of the ocean surface wave field. Details on this new parameterization, model initialization, and spinup, are described in Roach et al. (2019). Recalling the finite measurement resolution of the CryoSat2 dataset, the modeled representative radius is calculated only including floe size categories from 150 m (as floe sizes are radii, this corresponds to a radius equal to the minimum chord length) and larger. We include all FCD measurements here (chord lengths above 300 m) to make the broadest comparison, but note that the potential underrepresentation of floes with diameters near the sampling resolution may lead to inaccurate values of $\stackrel{\mathrm{\u203e}}{r}$ in regions mainly consisting of such floes.
Figure 6a–b, d–e compare modeled and observed climatologies of Arctic representative radius (for floe diameters 300 m and larger) averaged over 2011–2016 and the months of October–December (a, b) and February–May (c, d). Geographic variability of representative radius is broadly similar between model and observation: the largest floes lie in the Arctic interior, with regions of smaller floes in the straits and continental margins. Across the interior Arctic, simulated representative radii are significantly larger than are found in the observations, as the Roach et al. (2018a) FSTD model does not include processes that break up large floes in the absence of ocean surface waves. To compare seasonality between model and observations, we compare only those regions that experience wave fracture in the model runs, areas we collectively term the marginal ice zone (MIZ). The MIZ is defined by excluding categories that do not experience wave fracture in a given month (see Appendix D), shown as the contoured regions in Fig. 6a–b, d–f and greyed out in Fig. 6b, e. All excluded “pack ice” regions have modeled representative radii greater than 18 km. The MIZ region accounts for 37 % of grid areas with at least 25 chord measurements in months from October to December and 35 % of such areas for the period February–March. Note that the month of October is absent from these plots because no wellsampled regions are classified as MIZ across all model years according to the criteria outlined in Appendix D.
Figure 6c compares the observed (red) and modeled (blue) Arcticaverage representative radii for the MIZ over the period 2011–2016 as in Fig. 3a. The seasonal cycle of representative radius in the MIZ is different in the observations (red line; thin orange lines are individual months) than when all geographic regions are included (Fig. 3a). The seasonal cycle of representative radius in the pack ice region (i.e., not the MIZ) is shown as a green line in Fig. 6c. In the MIZ, average representative radii are smaller (on average 4.17 km vs. 6.49 km in the pack ice region). In contrast to the seasonal variation across all geographic regions (Fig. 3a) as well as in the pack ice, floes are larger in February–April than in November–December (5.40 km vs. 3.15 km). In both the MIZ and pack ice regions, however, average representative radius is similar in late winter. The largest difference between the two regions is from November to December, where representative radii are more than twice as large in the pack ice than the MIZ.
Figure 6f shows the annual average representative radius in the MIZ (red), pack ice (green), and modeled MIZ regions (blue). Modeled MIZ representative radii have a similar magnitude compared to the MIZ observations, though these regions have smaller floes than the interior. To address the scale mismatch between the toohigh modeled floe sizes and observed representative radii in the interior Arctic, as well as the strong and different seasonal cycle in representative radius in both regions, modeling efforts must include additional mechanisms for reducing floe size in the Arctic interior away from waves, such as mechanical fragmentation (Toyota et al., 2006; Rynders et al., 2016) or ridge dynamics (Roberts et al., 2019), to obtain realistic representative radii across the entire Arctic, as these processes are not present in the model used to make this comparison.
Here we developed and demonstrated a method for deriving the statistics of the sea ice FSD from satellite radar altimeter measurements of chord length. This method provides the first panArctic accounting of climaterelevant quantities derived from the FSD, permits testing of existing scaling laws previously used to characterize distributions of floe size, and allows for gridded comparisons between FSD models and observations. Using this new technique we produced climatological, annualaverage, and geographic mean moments of the Arctic FSD across a range of resolved length scales from 300 m to 100 km.
With the combination of satellite altimetry and mathematical theory, we were able to rigorously examine the powerlaw hypothesis related to the FSD under simple assumptions about the underlying floe chord data and the fidelity of CryoSat2 satellite retrievals. Segmenting measurements by geographic location, by month, and by year, we find limited statistical basis for a powerlaw scaling beginning below about 6.5 km. In a limited number of geographic locations, we find the observational data cannot rule out powerlaw scaling, except for typical sizes above about 6.5 km. Assuming a powerlaw floe size distribution can bias sea ice model output and conceptual understanding, the geographic variability and lack of consistent multiscale behavior reinforces the need for sea ice models to account for floescale processes rather than diagnose a distributional shape.
Observations that span the polar regions and different years and seasons are valuable for future refinement of processbased models of the FSD. In Sect. 5, we demonstrated how such model–observation comparisons can be made and can provide useful insights for model developers. At present, some general features of floe size evolution (in particular the magnitude and seasonal cycle of the representative radius) are broadly similar between model and observation in the marginal ice zone. Yet there is a significant scale mismatch in the interior Arctic between the presented simulations and this observational product because of missing fragmentation physics in the absence of ocean surface waves. Floe size modeling efforts have focused on marginal ice zone processes (Horvat and Tziperman, 2015; Zhang et al., 2015), and particularly floe sizes below about 1 km because these small floes play an important role in sea ice thermodynamics for floe sizes. The CryoSat2 observations, however, are best suited to resolving floe chords of several hundred meters and above. New satellite altimeters like ICESat2 have the potential to increase the chord length resolution to scales of 20–100 m and provide insight at smaller scales.
We emphasize strongly that refinement may be necessary to apply this method for operational purposes, trend analysis, and further model validation. This paper has focused on the framework for making altimetric measurements of the FSD and comparison to model output, but the obtained chord lengths and distributions have not been carefully validated against other observational methods, and this will be necessary before further application of this method. Before doing so, we have tried to outline the most significant uncertainties in the method. The typical assumptions of homogeneity, isotropy, and stationarity are invoked here at the length scale of the CICE model grid (O(25) km on each side) and timescale of 1 month. These statistical assumptions may not be satisfied if, for example, the number of measurements in a given region in 1 month is insufficient to sample the known anisotropy of the sea ice floe field, and additional passes change the mean chord length significantly (see Sect. S2 and Fig. S1). While we found little evidence for powerlaw scaling throughout most areas of the Arctic, this may be sensitive to the geographic (here the CICE model grid of approximately 25 km × 25 km) and temporal (here all measurements from 2010 to 2018) windows we use to collect and evaluate chord length measurements for a power law. The assumption of scaleinvariant sampling, observational uncertainty because of the finite sampling resolution, analysis of ambiguous returns, and the accuracy of retrievals in regions of thin sea ice may also affect the inferred size of sea ice floes. This in turn may affect the climatologies described in this study.
While processed CryoSat2 data have been validated against both visual imagery and groundbased observations, they were not designed with this application in mind – additional quality control may be necessary for climate studies of changing floe properties. The positive comparison between model and observation in Sect. 5 could also be due to a compensation between these measurement uncertainties and will need to be reexamined in future validation work. Yet observational uncertainties regarding, for example, the floe shape distribution can be roughly estimated at the order of the error in effective radius obtained for circular floes ($r=\sqrt{A/\mathit{\pi}}$) or a square ($r=\sqrt{A/\mathrm{4}}$), with a relative error of 25 %. Constraining model results beyond this scale of error will require further refinement. However, as shown in Fig. 6, at present the model–data mismatch in the interior Arctic can exceed a factor of 3. Even with expected levels of error in the present derived FCD–FSD product, some constraints on the model can be considered at present with this method. A future comparison of results from the IceSat2 and CryoSat2 altimeters will provide insights into the relevance of measurement and statistical uncertainties, as will comparison of altimetrically derived floe chord measurements with visual imagery.
Even accounting for important caveats that arise from making satellite measurements, remotely sensing the sea ice FSD from altimeters at subdaily resolutions can provide a significant increase in data for comparison and analysis of new sea ice models that parameterize the FSD. Previously the difficulty of making measurements of the FSD at relevant spatial and temporal scales has inhibited the widespread adoption of such floesensitive sea ice models. Understanding sea ice variability at the floe scale is also an important aspect of sea ice forecasting, and the ability to remotely assess the sea ice FSD at nearreal time will allow for further improvement of operational forecasting networks.
CPOM sea ice data, including raw floe length data, are available through the CPOM data portal at http://www.cpom.ucl.ac.uk/csopr/seaice.html (last access: 1 November 2019). The processed FCD–FSD statistics are available at https://github.com/chhorvat/CRYOSATFLOES/ (last access: 1 November 2019). The Roach et al. (2018a) FSTD model is publicly developed and available at https://github.com/lettieroach/ (last access: 1 November 2019).
For generic probability distributions S(D) and P(r), and a probability function $\stackrel{\mathrm{\u0303}}{F}(D;r)$, via Eq. (4) we have the relationship
where we restrict the upper bounds on the second integral because $\stackrel{\mathrm{\u0303}}{F}(D;r)$ is zero for D>r. Under the scaleinvariant sampling assumption $\stackrel{\mathrm{\u0303}}{F}(D;r)\mathrm{d}D=G\left(\mathit{\xi}\right)\mathrm{d}\mathit{\xi}$, where $\mathit{\xi}=\frac{D}{\mathrm{2}r}$ for D<2r (ξ<1). Therefore,
where A_{n} is the nth moment of G(ξ), a constant that depends on the functional form of G. For any such probability function (for example that derived in Sect. 2 for circular floes), the moments of the FSD and the moments of the FCD are proportional. Most of the hypothetical statistical distributions we would consider (for example, power laws) can be fully determined in terms of their moments, and thus the relationship between moments of the FSD and FCD is typically sufficient to reconstruct the underlying FSD.
Supposing P(r) were a powerlaw function, converting Eq. () to an integral over ξ from 0 to 1, we have
For a powerlaw function, $P(D/(\mathrm{2}\mathit{\xi}\left)\right)\propto {\left(\frac{D}{\mathit{\xi}}\right)}^{\mathit{\alpha}}$ and
From Eqs. (A4) and (A6), and under the assumptions of Sect. 2, all moments of the FSD and FCD are related by a computable function of the moment only, and powerlaw FSDs are derived from powerlaw FCDs with the same scaling law. While the proportionality of moments and Eq. (A6) prove that an observed powerlaw FCD must reflect an underlying powerlaw FSD, the same analysis used to arrive at Eq. (A6) can be repeated to find P(r) given a powerlawdistributed S(D) as well.
The real altimetric data product has a finite sampling resolution D_{min}, which can bias the computed FSD moments and powerlaw decay profile. For example, applied to real data with a finite sampling resolution, the integrals in Eqs. (4) to (5) are taken beginning at the minimum observed chord lengths D_{min} and floe sizes ${r}_{\mathrm{min}}={D}_{\mathrm{min}}/\mathrm{2}$. Moments of the distributions S and P reflect only statistics for floes larger than D_{min} and r_{min}, respectively. All other aspects of this derivation remain the same, as $\stackrel{\mathrm{\u0303}}{F}(D;r)$ is zero for any $r<D/\mathrm{2}$. However, the relationship expressed in Eq. (4) becomes
where $Y\left(r\right)\equiv {\mathrm{sin}}^{\mathrm{1}}\left(\frac{{D}_{\mathrm{min}}}{\mathrm{2}r}\right)$, ${S}_{n}\left(y\right)={\int}_{\mathrm{0}}^{y}{\mathrm{sin}}^{n}\left(x\right)dx$, and E is the error in relating the nth moments of S(D) and P(r). Since P(r) is unknown, E cannot be computed a priori. The function S_{n}(Y(r)) expresses the percentage of chords formed from floes of size r that would be smaller than D_{min}, although it is not readily expressed as a function of n. The most pathological distribution is when P(r) is a delta function at r_{min}, $P\left(r\right)=\mathit{\delta}(r{r}_{\mathrm{min}})$, $Y\left({r}_{\mathrm{min}}\right)=\mathit{\pi}/\mathrm{2}$, and E=1 as no chord lengths would be measured.
We can compute the error function for any delta function distribution as
and the misfit is the proportion of the integral of sin ^{n}(x) between 0 and Y(r^{*}). Because sin (x) is monotonically increasing from x=0 to π∕2, the integral of S_{n} is bounded above:
and the misfit error is bounded above by
The reciprocal of B is equal to π at n=0 and decreases sublinearly, and so away from r_{min} the error term decays exponentially with n and is small even for nearly pathological distributions (for n=1, ${r}^{*}={D}_{\mathrm{min}}$, for example, $E\le \mathit{\pi}/\mathrm{24}\approx \mathrm{14}\phantom{\rule{0.125em}{0ex}}\mathit{\%}$. Knowing the distribution of errors behaves in this way allows us to establish upper bounds by integrating P as a sum of δ functions.
We note that increasing resolution of floe chords will result in tighter bounds on this error. When $Y(r{)}^{*}\le \mathrm{1}$, which occurs when ${r}^{*}\ge \frac{{D}_{\mathrm{min}}}{\mathrm{2}\mathrm{sin}\left(\mathrm{1}\right)}\approx \mathrm{0.59}{D}_{\mathrm{min}}$, we can exploit a tighter bound using the fact that sin ^{n}(x)≤x^{n},
Using the same example as above (n=1, ${r}^{*}={D}_{\mathrm{min}}$) bounds the error $E\le {\mathit{\pi}}^{\mathrm{2}}/\mathrm{144}\approx \mathrm{7}\phantom{\rule{0.125em}{0ex}}\mathit{\%}$. A realworld distribution of floe sizes must have a peak value above zero; thus by increasing the sampling resolution (say, for example, to near the size of pancakes, i.e., D_{min}≈20 m or less, approached by the ICESat2 altimeter), this bound takes over and errors are reduced substantially.
We can explicitly solve Eq. (B3) for distributions with powerlaw tails. These distributions are peaked at the minimum floe size, and so will have high moment error. For power laws with $\mathit{\alpha}=\mathrm{1},\phantom{\rule{0.33em}{0ex}}\mathrm{2},\phantom{\rule{0.33em}{0ex}}\mathrm{3}$, or −4, $E\left(P\right(r;\mathit{\alpha},{r}_{\mathrm{min}}),\mathrm{1})$ is 1, 4, 16, or 25 %. For n=2, $E\left(P\right(r;\mathit{\alpha},{r}_{\mathrm{min}}),\mathrm{2})$ is .003, .04, 2, or 9.6 %: the increase in error with decreasing α is because sharper powerlaw slopes concentrate most of the distribution towards the smallest scale.
Given a set of floe chords {D}_{i} and an estimate of the beginning of a powerlaw tail D^{*}, we would like to find the most likely powerlaw floe size distribution $P(r;\mathit{\alpha},{r}_{\mathrm{min}})$ that generated them. As discussed in Appendix A, moments of the FSD and FCD are related by a multiplicative factor, and the distributions themselves will share the same powerlaw exponent. Thus we may test the powerlaw hypothesis directly on the FCD S(D). The powerlaw hypothesis means that S(D) is of the form
Following Muniruzzaman (1957) and Clauset et al. (2009) (see also the derivation in Stern et al., 2018a), we compute the loglikelihood of the observations for a given α (Eq. 10),
As the natural log is monotonically increasing in its argument, to find the most likely α, denoted $\widehat{\mathit{\alpha}}$, we take the derivative with respect to α and set to zero,
which resolves as a solution for the most likely α:
The above analysis concerns the most likely α that explains the FCD. We may ask a separate question: what is the most likely α, which we define as α_{P}, that would explain the FSD, given the explicit relationship that can be derived between S(D) and a powerlawdistributed P(r) examined in Eq. (10).
where C is unknown. Repeating the above analysis,
Next we take the derivative of ℒ with respect to α_{P} and setting to zero. We use the fact that $B(x,y)=B(y,x)$ and $\frac{\partial B(x,y)}{\partial x}=B(x,y)\left(\mathit{\psi}\left(x\right)\mathit{\psi}(x+y)\right)$, where ψ is the digamma function, to find
The maximum likelihood α_{P} is the solution to the transcendental equation,
which is an alternative method for obtaining the FSD scaling.
Due to limitations in the number of floe chords recorded at any particular location over time, we do not include all geographic locations when computing hemispheric means. Averaging is performed by including only geographic regions where there are at least 25 recorded floe chords. The area being averaged over is thus not fixed in time. For seasonal cycle plots, we only include months which have enough measurements for all fully sampled CryoSat2 years (2011–2018). For annual averages, we include only those years where all CryoSat2 months (excluding June–September) have enough measurements.
When masking additional regions to perform the model–observation comparisons in Fig. 6, we note that because the Roach et al. (2018a) model does not include processes that fragment larger floes into smaller floes in the absence of ocean surface waves, regions in the interior Arctic without wave activity have nearly all sea ice area belonging to the highest floe size categories. Nearly all regions where wave fracture is an active process also have representative radii below about 10 km (Roach et al., 2019). We define regions that do not experience wave fracture as those with an abnormally high simulated representative radius, which we choose to be the 22nd floe size category ($\stackrel{\mathrm{\u203e}}{r}=\mathrm{18.6}$ km) or above. The mask and comparisons in Fig. 6 are made by excluding all such areas.
The supplement related to this article is available online at: https://doi.org/10.5194/tc1328692019supplement.
CH derived the mathematical theory and wrote the paper. LR built and performed the climate model simulation. RT, AR, and AS provided and interpreted the CryoSat2 data. KH, CG, CB, and BK contributed to the study design. All authors have participated in paper preparation.
The authors declare no competing interests.
CH was supported by the NOAA Climate and Global Change Postdoctoral Fellowship Program, administered by UCAR's Cooperative Programs for the Advancement of Earth System Science (CPAESS), sponsored in part through cooperative agreement number NA16NWS4620043, years 2017–2021, with the National Oceanic and Atmospheric Administration (NOAA) and the U.S. Department of Commerce (DOC). CH, CG, and KH thank the American Mathematical Society for their support through the Mathematics Research Community “Differential Equations, Probability, and Sea Ice”, funded by NSF grants 1321794 and 1641020. LR was funded via Marsden contract VUW‐1408 and the New Zealand Deep South National Science Challenge, MBIE contract number C01X1445. CMB was supported by the National Science Foundation grant PLR1643431. BFK was supported by ONR grant N000141712963 and NSF grant 1350795. RT, AR, and AS were supported by the UK NERC Centre for Polar Observation and Modelling and the European Space Agency.
This research has been supported by the National Oceanic and Atmospheric Administration, Climate Program Office (grant no. NA16NWS4620043), the National Science Foundation, Division of Mathematical Sciences (grant nos. 1321794, 1641020, 1350795, and 1643431), the Office of Naval Research (grant no. N000141712963), the Marsden Fund (grant no. VUW‐1408), and the New Zealand Deep South National Challenge (grant no. C01X1445).
This paper was edited by Jennifer Hutchings and reviewed by Thomas Armitage and one anonymous referee.
Alberello, A., Onorato, M., Bennetts, L., Vichi, M., Eayrs, C., MacHutchon, K., and Toffoli, A.: Brief communication: Pancake ice floe size distribution during the winter expansion of the Antarctic marginal ice zone, The Cryosphere, 13, 41–48, https://doi.org/10.5194/tc13412019, 2019. a
Allard, R. A., Farrell, S. L., Hebert, D. A., Johnston, W. F., Li, L., Kurtz, N. T., Phelps, M. W., Posey, P. G., Tilling, R., Ridout, A., and Wallcraft, A. J.: Utilizing CryoSat2 sea ice thickness to initialize a coupled iceocean modeling system, J. Geophys. Res.Oceans, 62, 1265–1280, https://doi.org/10.1016/j.asr.2017.12.030, 2018. a
Armitage, T. W. K. and Davidson, M. W. J.: Using the Interferometric Capabilities of the ESA CryoSat2 Mission to Improve the Accuracy of Sea Ice Freeboard Retrievals, IEEE Trans. Geosci. Remote Sens., 52, 529–536, https://doi.org/10.1109/TGRS.2013.2242082, 2014. a
Armitage, T. W. K., Bacon, S., and Kwok, R.: Arctic Sea Level and Surface Circulation Response to the Arctic Oscillation, Geophys. Res. Lett., 45, 6576–6584, https://doi.org/10.1029/2018GL078386, 2018. a
Batchelor, G. K.: The theory of homogeneous turbulence, Cambridge University Press, 1953. a
Bateson, A. W., Feltham, D. L., Schröder, D., Hosekova, L., Ridley, J. K., and Aksenov, Y.: Impact of floe size distribution on seasonal fragmentation and melt of Arctic sea ice, The Cryosphere Discuss., https://doi.org/10.5194/tc201944, in review, 2019. a
Birnbaum, G. and Lüpkes, C.: A new parameterization of surface drag in the marginal sea ice zone, Tellus A, 54, 107–123, https://doi.org/10.1034/j.16000870.2002.00243.x, 2002. a
Clauset, A., Young, M., and Gleditsch, K. S.: On the Frequency of Severe Terrorist Events, J. Conf. Resolut., 51, 58–87, 2007. a, b
Clauset, A., Shalizi, C. R., Newman, M. E. J., Rohilla Shalizi, C., and J Newman, M. E.: PowerLaw Distributions in Empirical Data, SIAM Rev., 51, 661–703, https://doi.org/10.1137/070710111, 2009. a, b, c, d, e
Day, J. J., Hawkins, E., and Tietsche, S.: Will Arctic sea ice thickness initialization improve seasonal forecast skill?, Geophys. Res. Lett., 41, 7566–7575, https://doi.org/10.1002/2014GL061694, 2014. a
Feltham, D. L.: Sea Ice Rheology, Ann. Rev. Fluid Mechan., 40, 91–112, https://doi.org/10.1146/annurev.fluid.40.111406.102151, 2008. a
Gherardi, M. and Lagomarsino, M. C.: Characterizing the size and shape of sea ice floes, Sci. Rep., 5, 10226, https://doi.org/10.1038/srep10226, 2015. a
Guerreiro, K., Fleury, S., Zakharova, E., Kouraev, A., Rémy, F., and Maisongrande, P.: Comparison of CryoSat2 and ENVISAT radar freeboard over Arctic sea ice: toward an improved Envisat freeboard retrieval, The Cryosphere, 11, 2059–2073, https://doi.org/10.5194/tc1120592017, 2017. a
Herman, A.: Moleculardynamics simulation of clustering processes in seaice floes, Phys. Rev. E, 84, 1–11, https://doi.org/10.1103/PhysRevE.84.056104, 2011. a
Herman, A., Evers, K.U., and Reimer, N.: Floesize distributions in laboratory ice broken by waves, The Cryosphere, 12, 685–699, https://doi.org/10.5194/tc126852018, 2018 a
Horvat, C. and Tziperman, E.: A prognostic model of the seaice floe size and thickness distribution, The Cryosphere, 9, 2119–2134, https://doi.org/10.5194/tc921192015, 2015. a, b, c, d, e
Horvat, C. and Tziperman, E.: The evolution of scaling laws in the sea ice floe size distribution, J. Geophys. Res.Oceans, 122, 7630–7650, https://doi.org/10.1002/2016JC012573, 2017. a, b, c, d, e
Horvat, C. and Tziperman, E.: Understanding Melting due to Ocean Eddy Heat Fluxes at the Edge of SeaIce Floes, Geophys. Res. Lett., 45, 9721–9730, https://doi.org/10.1029/2018GL079363, 2018. a
Horvat, C., Tziperman, E., and Campin, J.M.: Interaction of sea ice floe size, ocean eddies, and sea ice melting, Geophys. Res. Lett., 43, 8083–8090, https://doi.org/10.1002/2016GL069742, 2016. a
Hunke, E. C., Lipscomb, W. H., Turner, A. K., Jeffery, N., and Elliott, S.: CICE : the Los Alamos Sea Ice Model Documentation and Software User's Manual Version 5.1 LACC06012, Tech. rep., Los Alamos National Laboratory, available at: http://oceans11.lanl.gov/trac/CICE (last access: 1 November 2019), 2015. a
Hwang, B., Wilkinson, J., Maksym, E., Graber, H. C., Schweiger, A., Horvat, C., Perovich, D. K., Arntsen, A. E., Stanton, T. P., Ren, J., Wadhams, P., Maksym, T., Graber, H. C., Schweiger, A., Horvat, C., Perovich, D. K., Arntsen, A. E., Timothy, P., Ren, J., and Wadhams, P.: Wintertosummer transition of Arctic sea ice breakup and floe size distribution in the Beaufort Sea, Elementa: Sci. Anthrop., 5, 40, https://doi.org/10.1525/elementa.232, 2017. a, b
JRA55: JRA55: Japanese 55year Reanalysis, Daily 3Hourly and 6Hourly Data, Research Data Archive at the National Center for Atmospheric Research, Computational and Information Systems Laboratory, 2013. a
Key, J.: Estimating the area fraction of geophysical fields from measurements along a transect, IEEE Tran. Geosci. Remote Sens., 31, 1099–1102, https://doi.org/10.1109/36.263782, 1993. a
Key, J. and Peckham, S.: Probable errors in width distributions of sea ice leads measured along a transect, J. Geophys. Res., 96, 18417, https://doi.org/10.1029/91JC01843, 1991. a
Kwok, R.: Arctic sea ice thickness, volume, and multiyear ice coverage: losses and coupled variability (1958–2018), Environ. Res. Lett., 13, 105005, https://doi.org/10.1088/17489326/aae3ec, 2018. a
Laxon, S., Peacock, N., and Smith, D.: High interannual variability of sea ice thickness in the Arctic region, Nature, 425, 947–950, https://doi.org/10.1038/nature02050, 2003. a
Laxon, S. W., Giles, K. A., Ridout, A. L., Wingham, D. J., Willatt, R., Cullen, R., Kwok, R., Schweiger, A., Zhang, J., Haas, C., Hendricks, S., Krishfield, R., Kurtz, N., Farrell, S., and Davidson, M.: CryoSat2 estimates of Arctic sea ice thickness and volume, Geophys. Res. Lett., 40, 732–737, https://doi.org/10.1002/grl.50193, 2013. a
Lindsay, R. W. and Rothrock, D. A.: Arctic sea ice leads from advanced very high resolution radiometer images, J. Geophys. Res., 100, 4533, https://doi.org/10.1029/94JC02393, 1995. a
Lüpkes, C. and Birnbaum, G.: Surface drag in the Arctic marginal seaice zone: A comparison of different parameterisation concepts, Bound.Lay. Meteorol., 117, 179–211, https://doi.org/10.1007/s1054600514458, 2005. a
Mandelbrot, B. B. and Wheeler, J. A.: The Fractal Geometry of Nature, Vol. 51, W. H. Freeman, https://doi.org/10.1119/1.13295, 1983. a
Manucharyan, G. E. and Thompson, A. F.: Submesoscale Sea IceOcean Interactions in Marginal Ice Zones, J. Geophys. Res.Oceans, 122, 9455–9475, https://doi.org/10.1002/2017JC012895, 2017. a
Muniruzzaman, A. N. M.: On Measures of Location and Dispersion and Tests of Hypotheses in a Pareto Population, Calc. Stat. Assoc. Bull., 7, 115–123, https://doi.org/10.1177/0008068319570303, 1957. a, b
Nere, N. K., Ramkrishna, D., Parker, B. E., Bell, W. V., and Mohan, P.: Transformation of the ChordLength Distributions to Size Distributions for Nonspherical Particles with Orientation Bias †, Indust. Eng. Chem. Res., 46, 3041–3047, https://doi.org/10.1021/ie0609463, 2007. a
Peacock, N. R. and Laxon, S. W.: Sea surface height determination in the Arctic Ocean from ERS altimetry, J. Geophys. Res.Oceans, 109, C07001, https://doi.org/10.1029/2001JC001026, 2004. a, b
Perovich, D. K. and Jones, K. F.: The seasonal evolution of sea ice floe size distribution, J. Geophys. Res.Oceans, 119, 8767–8777, https://doi.org/10.1002/2014JC010136, 2014. a, b
Pons, M.N., Milferstedt, K., and Morgenroth, E.: Modeling of chord length distributions, Chem. Eng. Sci., 61, 3962–3973, https://doi.org/10.1016/j.ces.2006.01.036, 2006. a
Quartly, G. D., Rinne, E., Passaro, M., Andersen, O. B., Dinardo, S., Fleury, S., Guillot, A., Hendricks, S., Kurekin, A. A., Müller, F. L., Ricker, R., Skourup, H., and Tsamados, M.: Retrieving Sea Level and Freeboard in the Arctic: A Review of Current Radar Altimetry Methodologies and Future Perspectives, Remote Sens., 11, 881, https://doi.org/10.3390/rs11070881, 2019. a
RichterMenge, J. A.: Relating arctic pack ice stress and deformation under winter conditions, J. Geophys. Res., 107, 8040, https://doi.org/10.1029/2000JC000477, 2002. a
Roach, L. A., Horvat, C., Dean, S. M., and Bitz, C. M.: An emergent sea ice floe size distribution in a global coupled ocean–sea ice model, J. Geophys. Res.Oceans, 123, 4322–4337, https://doi.org/10.1029/2017JC013692, 2018a. a, b, c, d, e, f, g, h, i
Roach, L. A., Smith, M. M., and Dean, S. M.: Quantifying growth of pancake sea ice floes using images from drifting buoys, J. Geophys. Res.Oceans, 123, 2851–2866, https://doi.org/10.1002/2017JC013693, 2018b. a
Roach, L., Bitz, C., Horvat, C., and Dean, S.: Advances in modelling interactions between sea ice and ocean surface waves, Journal of Advances in Modeling Earth Systems, 1–14, in review, 2019. a, b
Roberts, A. F., Hunke, E. C., Kamal, S. M., Lipscomb, W. H., Horvat, C., and Maslowski, W.: A Variational Method for Sea Ice Ridging in Earth System Models, J. Adv. Model. Earth Syst., p. 2018MS001395, https://doi.org/10.1029/2018MS001395, 2019. a
Rothrock, D. A. and Thorndike, A. S.: Measuring the Sea Ice Floe Size Distribution, J. Geophys. Res., 89, 6477–6486, https://doi.org/10.1029/JC089iC04p06477, 1984. a, b, c, d, e
Rynders, S., Aksenov, Y., Feltham, D., Nurser, G., and Naveira Garabato, A.: Modelling MIZ dynamics in a global model, in: EGU General Assembly Conference Abstracts, p. 1004, 2016. a
Sandberg Sorensen, L., Simonsen, S., Langley, K., Gray, L., Helm, V., Nilsson, J., Stenseng, L., Skourup, H., Forsberg, R., and Davidson, M.: Validation of CryoSat2 SARIn Data over Austfonna Ice Cap Using Airborne Laser Scanner Measurements, Remote Sens., 10, 1354, https://doi.org/10.3390/rs10091354, 2018. a
Schröder, D., Feltham, D. L., Tsamados, M., Ridout, A., and Tilling, R.: New insight from CryoSat2 sea ice thickness for sea ice modelling, The Cryosphere, 13, 125–139, https://doi.org/10.5194/tc131252019, 2019. a
Schulson, E. M. and Hibler, W. D.: The fracture of ice on scales large and small: Arctic leads and wing cracks, J. Glaciol., 37, 319–322, https://doi.org/10.1017/S0022143000005748, 1991. a
Skourup, H., Simonsen, S. B., Sorensen, L., Bella, A., Forsberg, R., Hvidegaard, S., and Helm, V.: ESA CryoVEx/EU ICEARC 2016Airborne Field Campaign with ASIRAS Radar and Laser Scanner over Austfonna, Fram Strait and the Wandel Sea, Tech. rep., National Space Institute, Technical University of Denmark, 2017. a
Smith, M. and Thomson, J.: Scaling observations of surface waves in the Beaufort Sea, Elementa: Science of the Anthropocene, 4, 000097, https://doi.org/10.12952/journal.elementa.000097, 2016. a
Squire, V. A.: Of ocean waves and seaice revisited, Cold Reg. Sci. Technol., 49, 110–133, https://doi.org/10.1016/j.coldregions.2007.04.007, 2007. a
Squire, V. A., Dugan, J. P., Wadhams, P., Rottier, P. J., and Liu, A. K.: Of Ocean Waves and Sea Ice, Annu. Rev. Fluid Mechan., 27, 115–168, https://doi.org/10.1146/annurev.fl.27.010195.000555, 1995. a
Steele, M.: Sea ice melting and floe geometry in a simple iceocean model, J. Geophys. Res.Oceans, 97, 17729–17738, https://doi.org/10.1029/92JC01755, 1992. a, b
Steer, A., Worby, A., and Heil, P.: Observed changes in seaice floe size distribution during early summer in the western Weddell Sea, Deep Sea Res. Part II, 55, 933–942, https://doi.org/10.1016/j.dsr2.2007.12.016, 2008. a
Stern, H. L., Schweiger, A. J., Stark, M., Zhang, J., Steele, M., and Hwang, B.: Seasonal evolution of the seaice floe size distribution in the Beaufort and Chukchi seas, Elem. Sci. Anth., 6, 48, https://doi.org/10.1525/elementa.304, 2018a. a, b, c
Stern, H. L., Schweiger, A. J., Zhang, J., and Steele, M.: On reconciling disparate studies of the seaice floe size distribution, Elem. Sci. Anth., 6, 49, https://doi.org/10.1525/elementa.304, 2018b. a, b, c, d, e, f
Tilling, R., Ridout, A., and Shepherd, A.: Nearrealtime Arctic sea ice thickness and volume from CryoSat2, The Cryosphere, 10, 2003–2012, https://doi.org/10.5194/tc1020032016, 2016. a, b
Tilling, R., Ridout, A., and Shepherd, A.: Estimating Arctic sea ice thickness and volume using CryoSat2 radar altimeter data, Adv. Space Res., https://doi.org/10.1016/j.asr.2017.10.051, 2018. a, b, c, d
Tilling, R., Ridout, A., and Shepherd, A.: Assessing the impact of lead and floe sampling on Arctic sea ice thickness estimates from Envisat and CryoSat‐2, J. Geophys. Res.Oceans, 124, https://doi.org/10.1029/2019JC015232, 2019a. a
Tilling, R., Ridout, A., and Shepherd, A.: Assessing the impact of lead and floe sampling on Arctic sea ice thickness estimates from Envisat and CryoSat‐2, J. Geophys. Res.Oceans, 124, https://doi.org/10.1029/2019JC015232, 2019b. a, b, c, d, e
Tolman, H. L. H. L.: User manual and system documentation of WAVEWATCH III TM version 3.14, Technical Note 276, 276, 194, available at: ftp://ftp.ifremer.fr/ifremer/cersat/products/gridded/wavewatch3/HINDCAST/publications/Tolman_etal_MMAB276_2009.pdf\%0Apapers3://publication/uuid/E1C39B58DBCB4F8AADCDF2D2210DDC46 (last access: 1 November 2019), 2009. a
Toyota, T., Takatsuji, S., and Nakayama, M.: Characteristics of sea ice floe size distribution in the seasonal ice zone, Geophys. Res. Lett., 33, L02616, https://doi.org/10.1029/2005GL024556, 2006. a, b
Toyota, T., Haas, C., and Tamura, T.: Size distribution and shape properties of relatively small seaice floes in the Antarctic marginal ice zone in late winter, DeepSea Res. Part II, 58, 1182–1193, https://doi.org/10.1016/j.dsr2.2010.10.034, 2011. a, b, c
Toyota, T., Kohout, A., and Fraser, A. D.: Formation processes of sea ice floe size distribution in the interior pack and its relationship to the marginal ice zone off East Antarctica, Deep Sea Res. Pt II, 131, 28–40, https://doi.org/10.1016/j.dsr2.2015.10.003, 2016. a
Tsamados, M., Feltham, D. L., Schroeder, D., Flocco, D., Farrell, S. L., Kurtz, N., Laxon, S. W., and Bacon, S.: Impact of Variable Atmospheric and Oceanic Form Drag on Simulations of Arctic Sea Ice*, J. Phys. Oceanogr., 44, 1329–1353, https://doi.org/10.1175/JPOD130215.1, 2014. a
Virkar, Y. and Clauset, A.: Powerlaw distributions in binned empirical data, The Ann. Appl. Stat., 8, 89–119, https://doi.org/10.1214/13AOAS710, 2014. a, b, c
Wadhams, P., Squire, V. a., Goodman, D. J., Cowan, A. M., and Moore, S. C.: The attenuation rates of ocean waves in the marginal ice zone, J. Geophys. Res., 93, 6799, https://doi.org/10.1029/JC093iC06p06799, 1988. a
Wernecke, A. and Kaleschke, L.: Lead detection in Arctic sea ice from CryoSat2: quality assessment, lead area fraction and width distribution, The Cryosphere, 9, 1955–1968, https://doi.org/10.5194/tc919552015, 2015. a
Wilchinsky, A. V. and Feltham, D. L.: Modelling the rheology of sea ice as a collection of diamondshaped floes, J. NonNewton. Fluid Mechan., 138, 22–32, https://doi.org/10.1016/j.jnnfm.2006.05.001, 2006. a
Wilchinsky, A. V. and Feltham, D. L.: Modeling Coulombic failure of sea ice with leads, J. Geophys. Res.Oceans, 116, C08040, https://doi.org/10.1029/2011JC007071, 2011. a
Williams, T. D., Bennetts, L. G., Squire, V. A., Dumont, D., and Bertino, L.: Waveice interactions in the marginal ice zone. Part 1: Theoretical foundations, Ocean Modell., 71, 81–91, https://doi.org/10.1016/j.ocemod.2013.05.010, 2013. a
Zhang, J., Schweiger, A., Steele, M., and Stern, H.: Sea ice floe size distribution in the marginal ice zone: Theory and numerical experiments, J. Geophys. Res.Oceans, 120, 3484–3498, https://doi.org/10.1002/2015JC010770, 2015. a, b
 Abstract
 Introduction
 Floe chords and the floe size distribution
 Climatology and trends in floe properties derived from CryoSat2 altimetry
 Evaluating the powerlaw hypothesis using floe size statistics derived from CryoSat2
 An example model–observation comparison of floe size variability
 Conclusions
 Data availability
 Appendix A: Proof that the FCD and FSD have the same statistical properties
 Appendix B: Bounds on the relationship between chord length and floe size moments
 Appendix C: Maximum likelihood estimation for chord length distributions
 Appendix D: Averaging and segmenting FSD statistics
 Author contributions
 Competing interests
 Acknowledgements
 Financial support
 Review statement
 References
 Supplement
 Abstract
 Introduction
 Floe chords and the floe size distribution
 Climatology and trends in floe properties derived from CryoSat2 altimetry
 Evaluating the powerlaw hypothesis using floe size statistics derived from CryoSat2
 An example model–observation comparison of floe size variability
 Conclusions
 Data availability
 Appendix A: Proof that the FCD and FSD have the same statistical properties
 Appendix B: Bounds on the relationship between chord length and floe size moments
 Appendix C: Maximum likelihood estimation for chord length distributions
 Appendix D: Averaging and segmenting FSD statistics
 Author contributions
 Competing interests
 Acknowledgements
 Financial support
 Review statement
 References
 Supplement