New methodology to estimate Arctic sea ice concentration from SMOS combining brightness temperature differences in a maximum-likelihood estimator

Monitoring sea ice concentration is required for operational and climate studies in the Arctic Sea. Technologies used so far for estimating sea ice concentration have some limitations, for instance the impact of the atmosphere, the physical temperature of ice, and the presence of snow and melting. In the last years, L-band radiometry has been successfully used to study some properties of sea ice, remarkably sea ice thickness. However, the potential of satellite Lband observations for obtaining sea ice concentration had not yet been explored. In this paper, we present preliminary evidence showing that data from the Soil Moisture Ocean Salinity (SMOS) mission can be used to estimate sea ice concentration. Our method, based on a maximum-likelihood estimator (MLE), exploits the marked difference in the radiative properties of sea ice and seawater. In addition, the brightness temperatures of 100 % sea ice and 100 % seawater, as well as their combined values (polarization and angular difference), have been shown to be very stable during winter and spring, so they are robust to variations in physical temperature and other geophysical parameters. Therefore, we can use just two sets of tie points, one for summer and another for winter, for calculating sea ice concentration, leading to a more robust estimate. After analysing the full year 2014 in the entire Arctic, we have found that the sea ice concentration obtained with our method is well determined as compared to the Ocean and Sea Ice Satellite Application Facility (OSI SAF) dataset. However, when thin sea ice is present (ice thickness . 0.6 m), the method underestimates the actual sea ice concentration. Our results open the way for a systematic exploitation of SMOS data for monitoring sea ice concentration, at least for specific seasons. Additionally, SMOS data can be synergistically combined with data from other sensors to monitor panArctic sea ice conditions.


Introduction
The Arctic Ocean is undergoing profound transformation.The rapid decline in Arctic sea ice extent and volume that is both observed and modelled (e.g.Comiso, 2012;Stroeve et al., 2012) may have become the key illustration of change in a warming planet, but change is widespread across the whole Arctic system (e.g.AMAP, 2012;IPCC, 2013;SEARCH, 2013).A retreating Arctic sea ice cover has a marked impact on regional and global climate, and vice versa, through a large number of feedback mechanisms and interactions with the climate system (e.g.Holland and Bitz, 2003;Cohen et al., 2014;Vihma, 2014).
The launch of the Soil Moisture and Ocean Salinity (SMOS) satellite, in 2009, marked the dawn of a new type of space-based microwave imaging sensor.Originally conceived to map geophysical parameters of both hydrological and oceanographic interest (e.g.Martin-Neira et al., 2002;Mecklenburg et al., 2009), SMOS is also making serious inroads in the cryospheric sciences (e.g.Kaleschke et al., 2010Kaleschke et al., , 2012;;Huntemann et al., 2014).Developed by the European Space Agency (ESA), the SMOS's single payload, called Microwave Imaging Radiometer using Aperture Synthesis (MI-RAS), is an L-band (1.4 GHz, or 21 cm wavelength) passive interferometric radiometer that measures the electromagnetic radiation emitted by Earth's surface.The observed brightness temperature (T B ) can be related to moisture content over the soil and to salinity over the ocean surface (Kerr et al., 2010;Font et al., 2013), which can be used to infer sea ice thickness (Kaleschke et al., 2012) and snow thickness (Maaß, 2013;Maaß et al., 2015).
Sea ice concentration (SIC), defined as the fraction of ice relative to the total area at a given ocean location, is often used to determine other important climate variables such as ice extent and ice volume.SIC has been the target of satellitebased passive microwave sensors such as the Special Sensor Microwave/Imager (SSM/I and SSMIS) and the Advanced Microwave Scanning Radiometer (AMSR-E and AMSR-2) for more than 30 years.SIC can be estimated due to the fact that the brightness temperature of sea ice and seawater are quite distinct.There exist a variety of algorithms with which to retrieve SIC from T B observations tuned to those higherfrequency sensors, that is, frequencies between 6 and 89 GHz (e.g.Cavalieri et al., 1984;Comiso, 1986;Ramseier, 1991;Smith, 1996;Markus and Cavalieri, 2000;Kaleschke et al., 2001;Shokr et al., 2008).Those algorithms present different advantages and drawbacks depending on frequency, spatial resolution, atmospheric effects, physical temperature, and others.According to Ivanova et al. (2015), the first source of error in the computation of sea ice concentration is the sensitivity to changes in the physical temperature of sea ice, in particular for those algorithms that use measurements between 10 and 37 GHz.They identify atmospheric water vapour and cloud liquid water as the second source of error except for algorithms at 89 GHz, where it becomes the dominant error.Another problem faced by higher-frequency radiometers is that the SIC retrievals are affected by the thickness of snow cover, which is difficult to determine.Wilheit (1978) analysed the sensitivity of microwave emissivity of open seawater to a variety of geophysical variables such as atmospheric water vapour, sea surface temperature, wind speed, and salinity as a function of frequency (Fig. 1).The figure illustrates that L-band (1-2 GHz) observations are in a sweet spot, with the effect of all variables but salinity being minimal around the SMOS frequency.The same author also showed that the signature of multi-year (MY) and first-year (FY) ice overlap in the lower microwave frequencies, while this is not the case at higher frequencies.
Although some authors (e.g.Mills and Heygster, 2011a;Kaleschke et al., 2013) have recently explored the feasibility of SIC determination using an aircraft-mounted L-band radiometer, a method that extends satellite-based SIC retrievals down to L-band (i.e.SMOS) frequencies has been missing.We therefore set out to develop a new method, which we present here.
A significant difference between high-frequency and Lband microwave radiometry is that ice penetration at the L band is non-negligible (Heygster et al., 2014).In other words, ice is more transparent (i.e.optically thinner) at low than at  Wilheit, 1978, andUlaby andLong, 2014).high microwave frequencies.As a consequence, the brightness temperature measured by an L-band antenna is emitted not only by the topmost ice surface layer but also by a larger range of deeper layers within the ice.Thanks to that increased penetration in sea ice (about 60 cm depending on ice conditions), the SMOS L-band radiometer is also sensitive to ice thickness (Kaleschke et al., 2012;Huntemann et al., 2014).
We exploit some of the SMOS observational features in this study to develop a new method to estimate SIC.These include a combination of acquisition modes involving dual and full polarization; continuous multiangle viewing between nadir and 65 • ; a wide swath of about 1200 km; spatial resolution of 35-50 km; and 3-day revisit time at the Equator or, more frequently, at the poles.In particular, the multiangle viewing capability of SMOS is a noteworthy feature; it means that the same location on the Earth's surface can be observed quasi-simultaneously from a continuous range of angles of incidence as the satellite overpasses it.
The new method we present in this paper uses SMOS brightness temperature observations T B and a maximumlikelihood estimator (MLE) to obtain SIC maps in the Arctic Ocean.We describe SMOS data and a radiative transfer model for sea ice that allows us to compute its emissivity in Sects. 2 and 3, respectively.We then introduce the concept of tie points and their sensitivity to different geophysical parameters to help with SIC retrievals via algorithmic inversion of SMOS data in Sects. 4.1,4.2,4.3,and 4.4, and the MLE inversion algorithm in Sect.4.5.We then perform an accuracy assessment of SIC estimates using SMOS by comparing them to an independent SIC dataset in Sect. 5 and with a discussion and conclusions in Sects.6 and 7, respectively.

SMOS data from the Arctic Ocean
Since the launch of SMOS in 2009, ESA has been generating brightness temperature full-polarization data products from it.In this study, we focus on the official SMOS Level 1B (L1B) product version 504 data north of 60 • N from 2014 to estimate SIC.The L1B data contain the Fourier components of T B at the antenna reference frame (Deimos, 2010), from which one can obtain temporal snapshots of the spatial distribution of T B (i.e. an interferometric T B image) by performing an inverse Fourier transform.The T B data are georeferenced at an Equal-Area Scalable Earth (EASE) Northern Hemisphere grid (Brodzik and Knowles, 2002) of 25 km on the side.The radiometric accuracy of individual T B observations from SMOS is ∼ 2 K at boresight, and it increases on the extended alias-free field of view (Corbella et al., 2011).Proceeding from L1B data, though computationally more demanding than the more traditional L1C data products, has several benefits.For example, it allows one to change the antenna grid from the operational size of 128 × 128 pixels to 64 × 64 pixels.As shown by Talone et al. (2015), the smaller grid is optimal in that it helps mitigate some of the spatial correlations between measurements that are present in the larger grid.
We correct T B for a number of standard contributions such as geomagnetic and ionospheric rotation and atmospheric attenuation (Zine et al., 2008).The galactic reflection is not significant at high latitudes, and no correction was applied.We then filter out outliers (defined as those estimates that deviate by more than 3σ from the mean value, where σ is the radiometric accuracy at the given point in the antenna plane).We also filter out T B observations in regions of the field of view that are known to have low accuracy due to aliasing (Camps et al., 2005), Sun reflections, and Sun tails.
To lower the noise level, we averaged T B measurements from both ascending and descending orbits over periods of 3 days, which thus define the time resolution of our SIC maps.We also averaged acquisitions in incidence angle binnings of 2 • .Since some incidence angles could be missing due to the SMOS acquisition feature and interferences, we use a cubic polynomial fit to interpolate T B measurements to have the full range of incidence angles in each grid position.

OSI SAF and other sea ice data products
We use SIC maps from the database of the Ocean and Sea Ice Satellite Application Facility (OSI SAF product version OSI-401a) of the European Organisation for the Exploitation of Meteorological Satellites (EUMETSAT) for comparison with the products we are obtaining.
These are computed from brightness temperature observations from SSMIS at 19 and 37 GHz, are corrected for atmospheric effects using forecasts from the European Centre for Medium-Range Weather Forecasts (ECMWF), use monthly dynamic tie points, are available on a polar stereographic 10 km grid for both polar hemispheres, and include SIC uncertainty estimates (Tonboe et al., 2016).In this study, we used daily SIC maps in the Arctic Ocean from the OSI SAF Northern Hemisphere products of the year 2014.
We also used SIC estimates from ice charts generated from various sensors by the National Ice Center (Fetterer and Fowler, 2009) to identify regions of interests to compute the 100 % ice tie points.

Theoretical model of sea ice radiation at microwave wavelengths
Passive radiometers measure brightness temperature T B at the antenna frame with different incidence angle.T B can be expressed as where ϒ is the atmosphere transmissivity, T B SURF the radiation emitted by the surface, T B ATM_DN the downward-emitted atmospheric radiation that gets scattered by the terrain in the direction of the antenna, and T B ATM_UP the upward-emitted atmospheric radiation.
The surface emission is defined as where θ is the incidence angle relative to zenith angle, e s the surface emissivity, and T the physical temperature of the radiation-emitting body layer.Hereafter, we will use T B to refer to surface brightness temperature, for simplicity.
The emissivity e and reflectivity of a layer are related by e = (1 − ).The reflectivity (sometimes also called R) is the ratio between reflected and incident radiation at the media boundaries for each polarization.for horizontal H and vertical V polarizations can be calculated using Fresnel equations, which depend non-linearly on the dielectric constant (ε) and on the incident θ i and refracted θ t angles: The frequency-dependent dielectric constant of a medium is a complex number defined as ε(f ) = ε (f ) + iε (f ), where the real part ε is related to the electromagnetic energy that can be stored in the medium, the imaginary part ε is related to the energy dissipated within the medium, and f is frequency.Note that brightness temperature varies linearly with emissivity (Eq.2) and hence also with reflectivity.
To calculate the brightness temperature T B of sea ice, we will assume a sea ice model consisting of horizontal layers of three media -air, snow, and thick ice.We use the incoherent approach (i.e.conservation of energy, instead of wave field treatment in the coherent approach).Then a plane-parallel radiative transfer model (Eq.4) is used to propagate to the surface the reflectivity computed at and through the ice-snow and snow-air media boundaries, making a number of simplifying assumptions.Specifically, our model assumes (a) that the media are isothermal and (b) that the thickness of the ice layer is semi-infinite so that radiation from an underlying fourth layer (i.e.seawater) does not need to be considered.This approach is similar to that used by other authors (e.g.Mills and Heygster, 2011b;Maaß, 2013;Schwank et al., 2015).These assumptions are realistic for the emission of sea ice that is thicker than about 60 cm at the observing frequency of SMOS, as discussed in Sect. 1, since the underlying seawater then makes no contribution to the overall emissivity.
To further simplify our approach, we assume that the snow layer in the model consists of dry snow, which is typical of winter Arctic conditions.Dry snow can be considered a lossless medium at 1.4 GHz, due to the fact that the imaginary part of ε is very small compared with the real part, as stated in Schwank et al. (2015).That means that there is no attenuation in the snow layer, and therefore its attenuation coefficient, α snow , is considered zero.We make this simplifying assumption because water in a wet snow layer would cause attenuation and therefore increase the total emissivity, but it is rarely possible to obtain meaningful data on the amount of water in wet snow.However, dry snow still has an effect on the refracted angle according to Snell's law, and hence on the emissivity, which is computed via Eq.( 3).The permittivity of dry snow depends on snow density (Tiuri et al., 1984;Matzler, 1996), which depends on the snow temperature.For a snow density of ρ s = 300 g cm −3 , the dry-snow permittivity at the L band is ε snow = 1.53 following the equation described in Schwank et al. (2015).
We can now define the simplified brightness temperature that results from an infinite number of reflections between the three medias as (Ulaby et al., 1986) where as and si are the reflectivity at the air-snow and snow-ice boundaries, respectively, and T snow and T ice are the physical temperature in the snow and ice layers, respectively.The term τ is the attenuation factor and is defined as τ = 2dα sec θ , where d is the depth of the snow layer and α the attenuation constant.T sky is the temperature of the cosmic background.The dependence of T B on θ , p, and f is embedded in the expressions of and τ .
The attenuation constant α of the middle layer, in the case of a low-loss medium (ε /ε 1), can be expressed as where c is the speed of light.The skin depth is defined as δ s = 1/α (m) and characterizes how deep an electromagnetic wave can penetrate into a medium (e.g.Ulaby and Long, 2014).
To compute the complex dielectric constant of sea ice ε ice , which is needed to compute si , we use the classic empirical relationship by Vant et al. (1978).In this model, permittivity depends linearly on the ice brine volume V b as where V br = 10 V b , and the coefficients a i can be obtained by linear interpolation to 1.4 GHz of the laboratory values from microwave measurements at 1 and 2 GHz (refer to Vant et al., 1978, for coefficient values).
The sea ice brine volume V b can be computed using Cox and Weeks (1983) as follows: where ρ, S, and T are sea ice density, salinity, and temperature, respectively.The F functions are cubic polynomials derived empirically, namely where the values of the coefficient a ij were given in Leppäranta and Manninen (1998) for ice temperatures between −2 and 0 • C, and for lower temperatures in Cox and Weeks (1983); see also Thomas and Dieckmann (2003).
Figure 2 shows the dependence of brightness temperature at the L band with angle of incidence for seawater and sea ice, as well as that of ice overlaid by a dry snow layer (following Eq. 4), for nominal Arctic temperature and salinity values.Specifically, temperature and salinity values used were after Maaß (2013): for seawater −1.8 • C and 30 g kg −1 , respectively, and for sea ice −10 • C and 8 g kg −1 .Note that the T B of seawater is significantly less than that of ice, and that the latter is slightly less than that of snow over ice.Also note the non-linear dependence of T B on incidence angle, the difference between H-and V-polarized waves for all three cases, and the larger variation with incidence angle of H than V over ice and snow (e.g.Maaß et al., 2015).
We also calculate the theoretical emissivity e s of a fourlayer model using the Burke et al. (1979) equation.The additional layer in this model is the seawater under sea ice, and we use the dielectric constant of seawater from Klein and Swift (1977).This layer does not need to be considered for the case of (optically) thick ice, but it becomes "visible" for the case of (optically) thin ice (i.e.thicknesses ≤ 60 cm, depending on ice temperature and salinity).The expression of T B for a four-layer model is defined in Burke et al. (1979) as where T i is the temperature of each layer, its reflectivity, γ the absorption coefficient, and z the layer thickness.The net effect of reducing the sea ice thickness and starting to sense seawater is a decrease in surface emissivity, and hence of T B (as illustrated in Fig. 5), relative to emissivity of thick ice (Shokr and Sinha, 2015).

Definition of robust indices from brightness temperature
It is rarely possible to obtain the ancillary geophysical data such as sea ice temperature, salinity, and ice thickness that are required to estimate brightness temperature from a microwave-emission model.Therefore, making assumptions and approximations becomes critically important.It is possible, however, to define a number of indices resulting from a combination of brightness temperature observations that are less sensitive to the unknown physical parameters.For example, estimates of soil moisture or sea ice concentration from radiometric measurements are often derived by combining T B measurements obtained from different polarizations, frequencies, and angles of incidence (Becker and Choudhury, 1988;Owe et al., 2001).
Hereafter, we use two indices, the polarization difference (PD) index and the angular difference (AD) index.The PD index is defined as the difference between T B measurements obtained at vertical T B V and horizontal T B H polarizations as The AD index is defined as the difference between two vertical polarization T B measurements obtained at two different angles of incidence as Figures 3 and 4 show the variation of PD and AD for the thick-ice model with angle of incidence, respectively.In defining AD, we use vertical rather than horizontal polarization because identification of the three media is facilitated by the larger dynamic range and non-crossing signatures of vertical polarization (Fig. 4).We choose θ = 35 • angle difference because this value represents a good compromise between sensitivity of the index and radiometric accuracy in the case of SMOS (Camps et al., 2005) and, importantly, is also well supported by the wide range of satellite viewing angles that characterizes SMOS.
Although the polarization ratio (PR) is also a commonly used index, we have chosen PD after verifying that its dynamic range is larger than that of PR and suspecting that PD would yield higher-accuracy estimates given the SMOS error budget.

Calibration of sea ice concentration using tie points
Tie points are widely used for retrieving SIC with higherfrequency radiometers, as well as in other fields such as pho-  The vertical line at 25 • incidence angle is drawn for reference to tie points, which are marked with a solid circle on vertical polarization for the three media.
togrammetry (e.g.Khoshelham, 2009).In this study, we use tie points as the typical T B values for 100 and 0 % concentrations, which permit us to compute the sea ice concentration.Tie points can therefore be viewed as SIC calibration points because their expected radiation can be unambiguously determined.
Figure 3 shows theoretical PD tie point values for open water and sea ice, as well as ice with a snow layer.The values for an angle of incidence of 50 • are marked by solid red circles.This angle represents a good compromise in PD contrast between the two media and SMOS accuracy (Camps et al., 2005).The two bounding values are 62.9 K for seawater and 26.8 K for ice with snow cover (Table 1).The large difference between tie point values suggests that it is possible to estimate SIC at the L band.
Figure 4 shows theoretical AD tie point values for difference in incidence angle θ = 35 • and angles of incidence up to θ = 30 • , which, per Eq. ( 11), represents the T B V difference between θ = 60 • and θ = 25 • .The values for an angle of incidence of 25 • are marked by solid red circles, for which the tie points are 51.8K for seawater and 8.6 K for ice with snow cover (Table 1).Hereafter, PD and AD are evaluated at the incidence angles of θ = 50 • and θ = 25 • , respectively.
Figure 5 shows that T B at nadir increases non-linearly as a function of ice thickness up to the saturation value of ∼ 250 K, which is reached when ice becomes about 70 cm thick.Notice that T B estimates start at an ice thickness of 5 cm because there is a discontinuity in the Burke model as the thickness of ice tends to zero (e.g.Kaleschke et al., 2010;Mills and Heygster, 2011a;Maaß, 2013;Kaleschke et al., 2013).Compared with T B , the total variation of both AD and PD with ice thickness is significantly smaller, and they are therefore better suited to estimate sea ice concentration.

Sensitivity of estimates of sea ice concentration to surface emissivity changes
In this section, we calculate the sensitivity of SIC estimates to changes in surface emissivity due to variations in the physical properties of sea ice (i.e.salinity, temperature, and thickness).We work with estimated SIC derived from the three indices T B , PD, and AD.This is done following a standard error propagation method (as also used in Comiso et al., 1997).
It is important to determine how changes in ice conditions affect SIC estimates through those three indices to try to minimize SIC errors obtained using SMOS.
Table 2 lists  It should be noted that those sensitivities are calculated using the model and the nominal Arctic temperature and salinity values defined in Sect.3. In order to assess which index is least sensitive to changes in a given geophysical variable, we calculate absolute sensitivities, defined as the sensitivities multiplied by the dynamic range of the measurements.
Knowing the value of the tie points of sea ice (SIC = 100 %) and seawater (SIC = 0 %), one can compute the average slopes of the SIC estimates to their corresponding parameters T B , PD, and AD (i.e.δSIC/δT B , δSIC/δPD, and δSIC/δAD).From data in Table 1, we obtain the average slopes as δSIC/δT B = 0.65, δSIC/δPD = 2.32, and δSIC/δAD = 2.77.These slopes can be used to propagate T B , AD, and PD errors to errors in the SIC estimates.
We assume reasonable values for the variability of the physical parameters on which our emissivity model depends -namely T , S, and d of ice (generically denoted by g) -as follows: T = 5 K, S = 4 g kg −1 , and d = 30 cm.Using the values in Table 2 and the calculated average slopes, one can finally compute the errors in SIC estimates associated with the geophysical variability of g when the index I is To evaluate the final impact of geophysical variability on the SIC evaluation using the index I , we compute the root sum square (RSS) of the SIC uncertainties due to the geophysical parameters (Table 3).The table shows that AD is the most robust index to retrieve SIC, slightly better than PD, and significantly better than T B .Because T B is theoretically more sensitive to thin ice than the other two indices, one can expect that the use of T B to retrieve SIC would result in larger SIC errors.Moreover, the uncertainty distribution of T B is too broad, especially due to thickness, and thus less adequate to fulfil the statistical hypotheses used to derive SIC.Despite the uncertainties in the theoretical physical model of ice, we consider the differences significant enough to focus on inversion algorithms using the PD and AD indices, and not on T B .

Comparison with empirical tie points
Following the theoretical analysis, we now turn to evaluating its performance empirically.We therefore select several regions of interest in the Arctic Ocean where SIC has been determined to be either 0 or 100 % by other sensors and methods.To identify such regions, we use SIC maps from OSI SAF and from the National Ice Center.In particular, we selected the open-seawater region within latitudes 55-70 • N and longitudes 20 • W and 25 • E, which comprises more than 2000 pixels in a typical SMOS image.For sea ice, we selected the MY ice region within latitudes 78-83 • N (the northernmost latitude observable by SMOS) and longitudes 75-150 • W, which comprises about 1000 pixels per SMOS image.We expect some level of uncertainty associated with the selection of the region to compute the 100 % tie points for summer periods stemming from known errors in the summer SIC products by OSI SAF (Tonboe et al., 2016).
We calculated SMOS brightness temperatures of these target regions to evaluate their potential as empirical tie points for seawater and sea ice.The temporal variation, in 2014, of the spatially averaged (median) T B at nadir of the two geographic regions is shown in Fig. 6  region, the figure shows that the brightness temperature is constant, at about 99 K, to within ∼ 2.5 K (1σ standard deviation) throughout the year.For the ice region, T B is also stable during the non-summer months, but it drops by about 20 K during the summer season due to changes in surface emissivity associated with snow and ice melt and concurrent formation of meltwater ponds.The factor-of-2 increase in formal error in summer relative to winter is also an indication of increased radiometric variability in surface conditions (as shown in Table 1).
Figure 7 shows that the temporal radiometric stability of the seawater region during 2014, and that of sea ice during the non-summer months, is also reflected in the AD and PD indices, as one would expect.This suggests that a different set of tie points during winter and summer periods could be beneficial for the quality of the SIC retrievals.On the other hand, the AD and PD tie point values are very stable during winter and spring (November to June), indicating that values are robust to variations in physical temperature and that it may not be necessary to compute tie point values often (daily or monthly), as done with the OSI SAF product.
Figure 8 shows a 2-D scatter plot of AD and PD indices for the two regions defined above during March (winter tie point) and July (summer tie point) 2014.The index values associated with seawater and with ice group form two welldifferentiated clusters, which implies that the two types of regions can be clearly segregated using these indices.This is also true for the summer tie points even though in this case the dispersion is larger and values are closer to sea tie points, as expected following Fig. 7.
The modelled (with snow and without) and observed TB, AD, and PD tie point values for winter and summer 2014, as well as the standard deviation (σ ) of the measurements, are listed in Table 1.It is encouraging that most of the values are in agreement at about 2σ despite underlying model assumptions such as uniform sea ice temperature and specular ocean surface.Another important result is that the observed SMOS data are closer to the model when snow is considered.

Retrieval of sea ice concentration
The brightness temperature of mixed pixels -that is, ocean pixels partially covered by sea ice -can be expressed as a linear combination of the brightness temperature of ice and seawater weighted by the percentage of each surface type (e.g.Comiso et al., 1997): where C is the fraction of ice present in a pixel, with C = 1 corresponding to 100 % of ice and C = 0 to 0 % of ice, or equivalently 100 % of seawater.Since AD and PD (Eqs.10-11) depend linearly on brightness temperature, Eq. ( 13) can be used to express both AD and PD.
There are several possible strategies to estimate sea ice concentration at a given pixel from the AD and PD values measured at that pixel.The simplest approach is to consider that the values of the tie points are good representatives of the values of AD and PD at the respective medium, i.e. seawater and sea ice, such that Concentration C can thus be retrieved from the value of either AD or PD by inverting the associated linear equation.In general, C can also be evaluated simultaneously with the AD and the PD observations by averaging the values obtained from both indices as This is known as the linear estimation of SIC.However, this approach might be too simple, as the values of AD and PD on ice and seawater can have some non-negligible dispersion due to geophysical conditions and to radiometric noise.
In this paper, a new inversion algorithm to estimate C is presented, which considers that AD and PD have known distributions, and by combining the observations it is possible to infer the value of C that is statistically more probable.
The distributions of the SMOS AD and PD are unimodal and symmetric (not shown), thus allowing us to approximate them by Gaussians and considering the pure ice and pure sea measurements as independent.Therefore we can easily use a MLE approach.The MLE has many optimal properties in statistical inference such as (e.g.Myung, 2003) sufficiency (the complete information about the parameter of interest is contained in the MLE estimator), consistency (the true value of the parameter that generated the data is recovered asymptotically, i.e. for sufficiently large samples), efficiency (asymptotically, it has the lowest-possible variance among all possible parameter estimates), and parameterization invariance (same MLE solution obtained independent of the parametrization used).
Assuming the linearity superposition of indices (Eq.14), it follows that the distributions of AD and PD (f AD , f PD ) in a general ocean pixel can be expressed as where the bar over the AD and PD indices refers to their mean values, the subindex identifies the medium, and σ is the associated standard deviation for each index and media.To obtain the mean and standard deviation values, we used the SMOS measurements at the regions for generating tie points and periods discussed in Sect.4.4.The symbol N means normal probability density function.
As a first approximation, we have considered AD and PD two independent variables.It thus follows that the likelihood function L is equal to the product of their distributions or, equivalently and conveniently, to their sum (recall that the likelihood is the logarithm of the probability density function): The MLE of SIC is the value of C that maximizes the likelihood function l.

Quality algorithm assessment
We have calculated AD and PD values from SMOS brightness temperature and used the MLE approach to obtain SIC estimates over the Arctic Ocean for the year 2014.We have estimated SIC using different tie points, characterized by their central value and dispersion.For seawater, we have used a single year-round median value and the associated standard deviation for each index.For ice tie points, we have used two sets of values, as suggested by the results in Fig. 7.For the first set, we have computed for all years the median of the tie points between December and May (Table 1), i.e. the winterspring months, when Arctic sea ice extent is close to its annual maximum.For the second set, we have used those same winter-spring values for the months of October through May but the average of the summer values for the months between June and September (Table 1).We have used neither the October nor the November data to compute ice tie point values because these are months of maximum extension of thin ice, and underlying emission through thin ice could cause some errors on the SIC estimates (Fig. 5 and Table 2).The root-mean-square (rms) error of SIC retrievals relative to OSI SAF over the Arctic Ocean is shown in Fig. 9. Four types of retrievals and two sets of tie points are compared.Introducing a specific set of summer tie points (black solid line) reduces the rms error with respect to using only one unique tie point for the whole year (black dotted line).The rms reduction is about 24 and 12 % in July and August, respectively, and to smaller degree in June and September.Therefore, we will hereafter use a different set of tie point values in summer and winter.Furthermore, using the set of summer-winter tie points, results from four types of inversions that stem from combinations of linear and MLE method and indices are compared in Fig. 9.The lowest rms values through all months in 2014 but January are obtained with the MLE inversion algorithm and the AD index alone.The evolution throughout the year of the rms obtained with the linear retrieval method is similar in the case of the MLE method, albeit with a ∼ 5-10 % increased noise level.Larger rms values and increased temporal variability are observed when the PD index is also used.The rms error of all retrievals is largest in fall, in particular if the PD index is used.Those are months of ice formation; therefore vast regions become covered with frazil ice, nilas, and thin young ice, following the minimum ice extension of September.All methods converge to similar results in September, since this period is the one with minimum ice extension and minimum thin ice is expected (thus resulting in very small difference between using AD or AD and PD methods).
The spatial variation of the difference in MLE SIC retrievals when using only the AD index and when using the AD and PD indices for the period 2-5 November 2014 is shown in Fig. 10.As expected, the largest differences are associated with regions of thin ice formation, in particular in the Laptev Sea, Kara Sea, and along the edge of the ice pack both in the western Arctic and in the Atlantic sector.Together, the spatio-temporal snapshots in Figs. 9 and 10 highlight the sensitivity of PD to the presence of thin ice, which naturally leads to an increase of the retrieval error when PD is used.This conclusion is not fully consistent with the analysis done using the models in Sect.4.3 on the dependence of the indices (T B , PD, AD) on ice thickness.Table 2 shows that, theoretically, PD is slightly less sensitive to thin ice than AD.However, the AD index is the least sensitive (lowest RSS) to variations of all the analysed variables.Therefore, we will

Accuracy assessment of SMOS SIC retrievals
We have evaluated the mutual consistency of the SMOS SIC retrievals, and in the process we have determined which is the approach that leads to the minimum error in the retrieval of SIC.We now evaluate the accuracy of those retrievals.Although a representative (in the space-time domain) groundtruth dataset that allows us to assess the accuracy of SMOS retrievals does not exist, the SIC estimates from OSI SAF are a good option for cross-check.They are independent from SMOS, the spatio-temporal sampling and resolution of their products are commensurate with those of SMOS, and their error budget is available.
The spatial distribution of SMOS SIC in the Arctic Ocean has been estimated from SMOS data for the 3-day period 2-5 March 2014, and it has been compared with the OSI SAF SIC product on 4 March 2014.The largest differences between both algorithms are located at the margins of the sea ice cover, where thinner ice can be expected (see Fig. 11).March is the month of maximum sea ice extent, but the results for other winter months are similar.
On the other hand, November is the month of maximum extension of thin young ice, especially through the Beaufort Sea since ice in the Laptev and Kara seas remains thin during winter (Shokr and Dabboor, 2013).Significantly larger differences between SMOS and OSI SAF products are now observed over a much wider area of the Arctic Ocean including the Barents, Kara, Laptev, East Siberian, and Beaufort seas (Fig. 12).
The brightness temperature measured by a passive microwave radiometer increases with sea ice thickness up to  a saturation value.Such an increase is more gradual for low frequencies and horizontal polarization (e.g.Ivanova et al., 2015).At the SMOS L band, the increase of emissivity with ice thickness reaches saturation for an ice thickness that is about 60 cm, depending on ice salinity and temperature (Kaleschke et al., 2012), whereas at the OSI SAF frequencies (19 and 37 GHz) it is only a few centimetres (Heygster et al., 2014;Ivanova et al., 2015).For example, for pixels that are 100 % covered by thin ice of, say, 25 cm thickness, the AD and PD values for those pixels will be slightly different than the tie point value of ice because the value of ice tie points was computed from thick MY ice (see Fig. 5) for model analysis.This contrast leads to a difference in classification of such pixels, which will be considered mixtures of water and ice in the case of SMOS and as 100 % ice with OSI SAF.In other words, the estimation of SIC of a sea covered by frazil ice and nilas will be higher for OSI SAF than for SMOS.
To further analyse this classification difference, we have calculated the probabilities of SMOS SIC conditioned by values of OSI SAF SIC using a full year, 2014, of Arcticwide estimates.The probability of estimating a SIC value with SMOS that is less than or equal to 5 % when the estimated OSI SAF SIC is 0 % is shown in Fig. 13 (red line).As expected, the conditioned probability is very high throughout the year.This implies that both products have a similar abil-  (1) OSI SAF SIC < 0.9, (2) OSI SAF SIC > 0.9 and SMOS SIC < 0.9, and (3) OSI SAF SIC > 0.9 and SMOS SIC > 0.9.
ity to detect (close to) 100 % ocean pixels.This implies that the probability of having high SMOS SIC values when OSI SAF is low is almost zero, which also means that the rate of triggering false alarms on ice detection with SMOS is low.However, the probability of estimating a SMOS SIC equal to or higher than 90 % while the OSI SAF SIC is 100 % is not constant during the year and decrease with respect to the previous case.During the winter period (between January and April), the conditioned probability is notably high (near 0.9) (see Fig. 13 blue line).Then it decreases sharply during spring and most notably in summer.This change in the conditioned probability starting in the spring could stem from a change in ice properties.As the snow becomes wetter with the onset of the melt season, the observed emissivity starts to change, and this varies with the observating frequency (different scattering response).The observed increase of the conditioned probability in June could be due to the use of a summer tie point (applied from June to September), which improved the rms with respect to OSI SAF as shown in Fig. 9.The low conditioned probability in fall can be explained by the presence of thin ice.
We have analysed the spatial distribution of the conditioned probability of SIC estimates for the months of March and November.Those regions where OSI SAF SIC is more than 0.9 while SMOS SIC is less than 0.9 (light blue colour in Fig. 14) outline the edge of the ice cover.This is in good correspondence with the expected areas of thin ice.Additionally, this condition is extended when analysing November data (Fig. 14b), when thin ice is more frequent in the Arctic.
During the winter months, the spatial coefficients of determination (r 2 ) between SMOS and OSI SAF SIC are high (more than 0.65), which again is consistent with our interpretation about the role of thin ice in SMOS SIC (see Fig. 15).As melt starts, the correlation between SIC estimates continues to be high, thanks to the use of the summer tie point in the retrieval.In September, ice cover extent is at minimum because ice growth has not started yet, there is almost no thin ice, and the correlation remains high.The correlation drops in the fall (between October and December) because ice growth starts by freezing of the sea surface, producing large amounts of new thin ice.To compute these values, we have only included SIC values between 0.05 (5 %) and 0.95 (95 %) when computing correlations to avoid the two extreme values leading to too-high, non-significant values of correlation.

Discussion
The two PD and AD indices, which are derived from brightness temperature, have been designed to maximize their differences between open water and sea ice.Both have a low response to changes in the geophysical characteristics of the media, which has been confirmed by using theoretical models and by performing sensitivity analysis.
The tie points, defined as the characteristic values of our reference indices on the different media, have been calculated from SMOS data.Compared to theoretical values, some small discrepancies at the 10-20 % level have been observed, probably due to simplifying assumptions such flat surface ice, flat sea, and constant temperature at the layers used in the theoretical models.We have thus decided to follow a more empirical approach.The use of two sets of tie points, one for summer and one for winter measurements, improves the results of the summer SIC maps relative to static unique tie points.This improvement is not caused by changes in the ice or sea physical temperature but most probably by changes in the ice properties, because as snow and ice become wetter during the melt season, the observed radiometric emission changes.This effect is also observed in measurements from radiometers at higher frequencies than SMOS.
We have introduced the MLE inversion algorithm to retrieve SIC from SMOS data.The method is based on the maximization of the a posteriori likelihood of the joint distribution of AD and PD indices, assuming that they are independent and normally distributed.This MLE algorithm is more robust (less noisy) than the linear inversion (Eq.14).It also improves the retrieved SMOS SIC with respect to a linear inversion method because the former takes into account the dispersion (error) of the tie points (reference), which makes the algorithm more robust to T B errors.SIC maps obtained using only the AD index are of better quality than when the AD and PD indices are used together.We attribute this to the higher sensitivity of PD than AD to physical changes in the media.
SMOS and OSI SAF SIC maps compare well in terms of correlation (determination coefficient higher than 0.65) and rms except in areas of thin sea ice.This difference can be explained by the higher penetration of SMOS in sea ice (about 60 cm) than that from higher-frequency radiometers.Thus, when ice is thinner than 60 cm, SMOS data lead to lower values of SIC that OSI SAF.

Conclusions
Estimating SIC using L-band observations such as those from SMOS is recommended for the negligible effect of the atmosphere on brightness temperatures, and also because the vertical polarization of T B is insensitive to snow depth (Maaß et al., 2015).
Two indices derived from brightness temperature, PD and AD, have been chosen, since they verify the two required conditions: they maximize the difference between open water and sea ice, and they present a low response to changes in the geophysical characteristics of the media.AD and PD tie point values have been shown to be very stable during winter and spring periods (Fig. 7), indicating that the values are robust to variations in physical temperature.Thanks to that, one can safely assume two sets of static (i.e.not temporally varying) tie points, one for summer and one for winter for SMOS data, and not fortnightly or monthly as is done in the case of the OSI SAF product.
We have shown that the best configuration for SIC retrieval is using AD only with the MLE inversion method.We exclude PD and T B because they are more sensitive to ice thickness; therefore, the combined use of AD and PD presents larger errors when thin ice is present (fall).The MLE inversion method presents better results than a linear inversion since it takes into account the uncertainty of the tie points.
SIC estimates from SMOS have some drawbacks with respect to those from higher-frequency radiometers.For example, whereas the spatial resolution of the high-frequency SIC estimates can reach ∼ 3 km, the resolution from SMOS will not be better than about 35 km.A second issue of SMOS is that it underestimates SIC in the presence of thin ice, which is characteristic of the ice edges and freeze-up periods.Therefore, SMOS data should be used in combination with some form of spatial masking for those regions.We suggest that SIC estimates from SMOS can complement those from higher-frequency radiometers, together yielding enhanced SIC products.
This dataset could be very beneficial during the summer period, since SMOS SIC, theoretically, should be less sensitive to summer metamorphosis, due to the larger wavelength.Previous works show that the T B and SIC measured at 6.9 GHz band are more robust to summer ice changes than higher-frequency measurements (Kern et al., 2016;Gabarro, 2017).The confirmation of this statement will be done as future work.
Data availability.Level 1 SMOS data can be obtained from the ESA distribution system (https://earth.esa.int/web/guest/data-access).OSI SAF sea ice concentration is provided by OSI SAF via ftp://osisaf.met.no/archive.Sea ice charts are generated by the National Ice Center (https://nsidc.org/data/g02172).
Author contributions.CG and PE conceived the study.CG and JAP-R (during his master's thesis) developed the study to choose the best indices using modelling simulations.PE gave advice on physics of the ice and on the datasets available for validation.AT contributed with new mathematical concepts, and MP gave support in the sensitivity analysis.CG wrote the paper, and all co-authors contributed to the discussion and gave input for writing.

Figure 1 .
Figure 1.Sensitivity of brightness temperature for open seawater over a range of observing frequencies in the microwave band for a set of key geophysical parameters (created after Wilheit, 1978, and Ulaby and Long, 2014).

Figure 2 .
Figure 2. Theoretical variation of brightness temperature with angle of incidence at the L band for (blue) seawater, (black) sea ice, and (red) a snow layer overlying a sea ice layer for (continuous) horizontal and (dashed) vertical polarizations.

Figure 3 .
Figure3.Modelled variation of polarization difference (PD) index with angle of incidence for (blue) seawater, (black) sea ice, and (red) a snow layer overlying a sea ice layer at the L band.The vertical line at 50 • incidence angle is drawn for reference to tie points, which are marked with a solid circle for the three media.

Figure 4 .
Figure 4. Modelled variation of angular difference (AD) index with angle of incidence for (blue) seawater, (black) sea ice, and (red) a snow layer overlying a sea ice layer for (continuous) horizontal and (dashed) vertical polarizations, and for θ = 35 • , at the L band.The vertical line at 25 • incidence angle is drawn for reference to tie points, which are marked with a solid circle on vertical polarization for the three media.

Figure 5 .
Figure 5. Theoretical variation with sea ice thickness of (blue; left axis) T B at nadir, (green; right axis) polarization difference (PD) at 50 • incidence angle, and (red; right axis) angular difference (AD) at θ = 25 • after the model by Burke et al. (1979) for a sea ice salinity of 8 g kg −1 , a sea ice temperature of −10 • C, and a snow layer of 10 cm thick over the ice.Note the factor-of-10 change between the left and right vertical scales.

Figure 6 .
Figure 6.Temporal variation of the average brightness temperature T B at nadir for (top) multi-year sea ice and (bottom) seawater at the two regions for generating tie points.Note the factor-of-4 change in the vertical scales.

Figure 7 .Figure 8 .
Figure 7. Same as Fig. 6 except here for angular and polarization difference indices.

Figure 9 .
Figure 9.Comparison against OSI SAF of one tie point (black dotted line) vs. two tie points (black plane line) with MLE, and MLE vs. linear retrieval techniques.If not defined in the labels, it is two tie points.

Figure 10 .
Figure 10.SMOS SIC with MLE AD+PD minus SMOS SIC with MLE AD inversion techniques for 3 November 2014.SIC scale is presented from 0 to 1.

Figure 13 .
Figure 13.Probability of having SMOS SIC more than 0.90 where OSI SAF SIC = 1 (blue line) and SMOS SIC less than 0.05 where OSI SAF SIC = 0 (red line) for 2014.Summer tie points are used for retrievals from June to September.

Figure 15 .
Figure 15.Coefficient of determination (R 2 ) between SMOS and OSI SAF SIC for 2014, considering only SIC data in the range from 5 to 95 %.

Table 1 .
Modelled (with and without snow) and SMOS-observed T B , PD, and AD median values.Errors quoted are the standard deviation around the median.

Table 2 .
Sensitivity of measurement T B , PD, and AD to ice temperature (T ), salinity (S), and thickness (d).

Table 3 .
Propagated SIC error using each index, computed from Eq. (12) for assumed (T , S, d) variations and root sum square (RSS).