Surface Properties Linked to Retrieval Uncertainty of Satellite Sea-Ice Thickness with Upward-Looking Sonar Measurements

One of the key sources of uncertainties in sea ice freeboard and thickness estimates derived from satellite radar altimetry results from changes in sea ice surface properties. In this study, we analyse this effect, comparing upward-looking sonar (ULS) measurements in the Beaufort Sea over the period 2003–2018 to sea ice draft derived from Envisat and Cryosat-2 data. We show that the sea ice draft growth underestimation observed for the most of winter seasons depends on the surface properties preconditioned by the melt intensity during the preceding summer. The comparison of sea ice draft time series in the Cryosat-2 era indicates that applying 50% retracker thresholds, used to produce the European Space Agency’s Climate Change Initiative (CCI) product, provide better agreement between satellite retrievals and ULS data than the 80% threshold that is closer to the expected physical waveform interpretation. Our results, therefore, indicate compensating error contributions in the full end-to-end sea-ice thickness processing chain, which prevents the quantification of individual factors with sea-ice thickness/draft validation data alone.


Introduction
Satellite radar altimeter measurements are widely used to retrieve sea ice thickness in the Arctic region at the basin-wide scale-e.g., [1][2][3][4][5]. The applied methods are based on sea ice freeboard retrieval under assumption that Ku-band radar signal backscatter from the snow-ice interface. The sea ice freeboard is defined as a height of the ice surface above the local sea level and determined from comparison of elevation measurements over sea ice floes and leads-i.e., cracks in the ice. The freeboard measurements are then converted to sea ice thickness, assuming hydrostatic equilibrium of floating ice. Therefore, uncertainties in sea ice thickness (SIT) estimates stem from freeboard retrieval and freeboard to thickness conversion. One of the main sources of uncertainties in sea ice freeboard retrieval is associated with snow properties [6][7][8]. Depending on snow temperature and salinity, the scattering horizon may migrate within the snow pack, which impacts the surface elevation estimates. Another uncertainty is related to the possible underestimation of the sea level, since specular returns identified as the measurements in leads dominate the radar echo even being hundreds meters off-nadir [9]. Increased surface roughness associated with ridges and hummocks can be the reason for one more error source as it impacts the scattering surface elevation and prevents the accurate identification of the mean sea ice surface [10]. The largest uncertainties in sea ice freeboard to thickness conversion are related to the lack of knowledge about snow depth and sea ice density. Snow climatology, based on measurements of drifting stations during the second half of 20th century [11], is still used as the main source of snow depth and density data. Variability of the sea ice density is also poorly known, and the fixed values evaluated-e.g., in [12]-are usually used for the sea ice thickness retrieval.
Thus, the validation of different retrieval algorithms is highly important for the evaluation of the sea ice thickness product quality, as well as the future evolution of sea-ice thickness algorithms. Sea ice draft measurements of the moored upward-looking sonars (ULS) [13][14][15] and Operation IceBridge (OIB) airborne measurements of sea ice freeboard and snow depth [16] are the main data sources used to validate satellite altimeter-based sea ice thickness. OIB airborne measurements were acquired every year in late spring by lidar and snow radar and can be used to perform direct validation of the sea ice freeboard retrievals along the airborne tracks over the western part of the Arctic. The ULS data from the Beaufort Sea Exploration Project (BGEP) provide continuous measurements of the sea ice draft-i.e., portion of the sea ice below the sea level, over mooring sites and allow the validation of temporal changes of sea ice thickness estimates. Since draft observations represent the main portion of the sea ice thickness, the ULS data allow an end-to-end validation of satellite SIT product, but do not provide information about sea ice freeboard.
In this work we use BGEP ULS measurements to analyse satellite radar altimeter retrievals derived by the algorithm described in [3] and applied to produce SIT data records in the frame of the European Space Agency's Sea Ice Climate Change Initiative (CCI). The CCI Phase 2 (CCI-2) SIT product was validated by Kern et al. [17], where the comparison of the gridded sea ice draft estimated using radar altimeter observations from Envisat and Cryosat-2 (CS2) satellites with, in particular, the BGEP ULS measurements, is provided. A strong correlation between CS2 retrievals and BGEP ULS measurements is reported in several studies-e.g., in [5,[18][19][20]. In continuation of this work, we present a detailed description of seasonal and inter-annual changes in the discrepancies between satellite-derived estimates and the BGEP ULS data. Further, the origin of these discrepancies is analysed in conjunction with the satellite and geophysical parameters.

Radar Altimeter Data
Sea ice freeboard and thickness are retrieved from the Envisat pulse-limited radar altimeter measurements (2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012) and SAR mode beam-limited Cryosat-2 (CS2) radar altimeter data (2010-2019) by the approach that minimizes seasonally dependent inter-mission bias [3,21,22]. Sea ice freeboard is estimated from the difference between elevations measured in leads and over sea ice surface by applying the threshold first-maximum retracker algorithm, which uses a percentage threshold of the waveform's first-maximum power return to identify retracking point [23]. While the fixed 50% threshold is used to retrack CS2 waveforms both for lead and sea ice surfaces, an adaptive threshold, which depends on waveform shape parameters, is applied to the Envisat measurements over sea ice to form consistent sea ice freeboard time series [3]. The choice of the threshold percentage is shown to affect the retrieved elevation, depending on surface properties [10,23]. Therefore, to evaluate the effect of the freeboard retrieval algorithm, we computed sea ice freeboard and corresponding thickness by applying an 80% threshold to the CS2 data over leads and sea ice floes. For the conversion of sea ice freeboard to thickness, we used snow climatology [11] modified over first-year ice [24] and ice-type dependent sea ice density [12] in accordance with the CCI algorithm [25].
Since wet snow and the formation of melt ponds in summer prevents the retrieval of sea ice freeboard by radar altimetry, the satellite estimates are limited to the period from October to April. For comparison with the ULS measurements, we calculated monthly mean sea ice freeboard and thickness using the along-track estimates within the 200 km radius of the mooring locations following the approach used in [18].

Upward-Looking Sonar Data
The Beaufort Gyre Exploration Project (BGEP) presently provides the sea ice draft  measurements  for  15  winter  seasons  from  2003-2004 to 2017-2018 (https://www.whoi.edu/beaufortgyre). The measurements from the four (A, B, C, and D) moored ULS located in the Beaufort Sea ( Figure 1) are available for different time periods (Table 1). Data  from the mooring A are the only ULS measurements available for the entire period of observations,  while there are data gaps during the winter seasons 2005-2006 and 2009-2010 for the mooring B, and  data from the moorings C and D are available to the season 2007-2008 and from the season 2006-2007, respectively. Thus, data from the moorings A, B, and D are available for the entire period of  CS2 measurements, while the moorings providing measurements over the period of Envisat  observations vary from year to year. The ULSs operate from around 50 m below the surface with a frequency 420 kHz (beamwidth 1.8 deg) and measures sea ice draft every 2 s within the 2 m footprint with the accuracy of individual range measurement of ± 5 cm. For comparison with the estimates derived from the radar altimeter data we calculated monthly mean sea ice draft, excluding measurements of less than 10 cm, as they may represent leads.
The ULSs operate from around 50 m below the surface with a frequency 420 kHz (beamwidth 1.8 deg) and measures sea ice draft every 2 s within the 2 m footprint with the accuracy of individual Remote Sens. 2020, 12, 3094 4 of 16 range measurement of ±5 cm. For comparison with the estimates derived from the radar altimeter data we calculated monthly mean sea ice draft, excluding measurements of less than 10 cm, as they may represent leads.

Satellite-Derived vs. ULS-Measured Sea Ice Draft
A monthly averaged sea ice draft, estimated as a difference between sea ice thickness and freeboard retrieved from satellite altimeter data, was compared with the BGEP ULS sea ice draft measurements. We consider two satellite sea ice draft monthly time series: (1) related to the mooring A, and (2) calculated as a mean of sea ice draft monthly values related to all moorings where ULS data are available during the specific winter season. The availability of the ULS measurements from the mooring A over almost the entire period of satellite observations gives us a consistent time series, while combining data from all moorings provides time series that are less noisy than those based on just one mooring. Figure 2 presents a comparison of the ULS measurements with three satellite-derived sea ice draft estimates: one from Envisat (dr Env ) and two from CS2 with elevations retrieved using 50% (dr CS2_50 ) and 80% (dr CS2_80 ) retracker thresholds. The differences between the ULS sea ice draft measurements (dr uls ) and consistent time series of dr Env and dr CS2_50 are comparatively small on average. The mean (root-mean-square (rms)) of sea ice draft differences, dr Env -dr uls and dr CS2_50 -dr uls , are −0.12 m (0.38 m) and 0.07 m (0.21 m) when comparing data for the mooring A, and −0.18 m (0.34 m) and 0.08 m (0.19 m) when combining data from all moorings. The 80% threshold on the leading-edge of the waveform corresponds to larger range and consequently lower freeboard estimate, since the higher threshold has a larger impact over sea ice waveforms than the narrow lead waveforms. Therefore, the estimates of sea ice draft dr CS2_80 are, as expected, lower than those of dr CS2_50 with the mean (rms) differences dr CS2_80

Satellite-Derived vs. ULS-Measured Sea Ice Draft
A monthly averaged sea ice draft, estimated as a difference between sea ice thickness and freeboard retrieved from satellite altimeter data, was compared with the BGEP ULS sea ice draft measurements. We consider two satellite sea ice draft monthly time series: (1) related to the mooring A, and (2) calculated as a mean of sea ice draft monthly values related to all moorings where ULS data are available during the specific winter season. The availability of the ULS measurements from the mooring A over almost the entire period of satellite observations gives us a consistent time series, while combining data from all moorings provides time series that are less noisy than those based on just one mooring. Figure 2 presents a comparison of the ULS measurements with three satellite-derived sea ice draft estimates: one from Envisat (drEnv) and two from CS2 with elevations retrieved using 50% (drCS2_50) and 80% (drCS2_80) retracker thresholds. The differences between the ULS sea ice draft measurements (druls) and consistent time series of drEnv and drCS2_50 are comparatively small on average. The mean (root-mean-square (rms)) of sea ice draft differences, drEnv-druls and drCS2_50-druls, are −0.12 m (0.38 m) and 0.07 m (0.21 m) when comparing data for the mooring A, and −0.18 m (0.34 m) and 0.08 m (0.19 m) when combining data from all moorings. The 80% threshold on the leading-edge of the waveform corresponds to larger range and consequently lower freeboard estimate, since the higher threshold has a larger impact over sea ice waveforms than the narrow lead waveforms. Therefore, the estimates of sea ice draft drCS2_80 are, as expected, lower than those of drCS2_50 with the mean (rms) differences drCS2_80-druls of −0.42 m (0.47 m) and −0.44 m (0.49 m).    satellite retrievals and ULS measurements are revealed from time series of the differences between these datasets (Figure 2c,d). In particular, the winter sea ice growth is underestimated for the most of seasons, as represented by the decrease in differences during winter. The mean seasonal change rate of the difference dr Env -dr uls is −7.2 cm/month and −4.9 cm/month for the mooring A and for combined data from all moorings, respectively. For the CS2 estimates, the larger underestimation of seasonal sea ice draft growth is observed for the dr CS2_80 draft: mean seasonal change rate of the difference dr CS2_80 -dr uls is −7.6 cm/month for the mooring A and −6.5 cm/month for all moorings. A corresponding seasonal change in the difference dr CS2_50 -dr uls is on average −4.8 cm/month and As well as the underestimation of seasonal growth the inter-annual variability of the sea ice draft differences observed during the CS2 observation period is larger for the dr CS2_80 estimates, as compared to the dr CS2_50 data. Moreover, the time series obtained using combined data from all moorings ( Figure 2d) suggest that there is a relationship between inter-annual changes and change rates of the sea ice draft difference, especially for the dr CS2_80 -dr uls values. In particular, the lowest seasonal mean sea ice draft differences during the seasons 2013-2014 and 2014-2015 correspond to the lowest underestimation of seasonal draft growth.

Waveform Parameters and Sea Ice Freeboard Retrievals
The observed seasonal and inter-annual discrepancies between satellite-derived and ULS-measured sea ice draft can result from the convoluted uncertainties in sea ice freeboard retrieval and freeboard to thickness conversion. Aiming to reveal the effect of the former, we formed monthly time series of the radar altimeter waveform parameters-leading edge width (LeW) and pulse peakiness (PP) derived, following [3,23], which can be used to characterize the sea ice surface properties in the vicinity of the moorings ( Figure 3). For the CS2 observation period, we also calculated the difference fb CS2_80 -fb CS2_50 -i.e., the difference between monthly mean sea ice freeboards used for estimating drafts dr CS2_50 and dr CS2_80 (Figure 2c,d). In the following sections, we use this difference as a measure for freeboard retrieval uncertainty.
For most of the winter seasons, an increase in the LeW and decrease in PP during winter reflect a change in the radar altimeter echo waveforms to a more diffuse shape. This represents primarily a sea ice surface roughness increase with an increase in the fraction of deformed ice and pressure ridges, and also an increase in volume scattering due to snow depth growth. Sea ice freeboard difference fb CS2_80 -fb CS2_50 is, as expected, always negative, as the lower height of sea ice floes are estimated with the 80% retracker threshold. During winter seasons the fb CS2_80 -fb CS2_50 values typically becomes more negative following changes in the LeW and PP parameters [23]. Seasonal changes in the differences fb CS2_80 -fb CS2_50 can be attributed to an increase in the distance between 50% and 80% thresholds for the waveforms with a more flattened leading edge. This is consistent with the conclusion reported by Landy et al. [10], in that location on the leading edge, corresponding to the mean surface elevation within the footprint, depends on the sea ice surface roughness. Seasonal decrease in fb CS2_80 -fb CS2_50 in the area of BGEP moorings location is not observed only during the 2013-2014 and 2014-1015 seasons, which is in agreement with small seasonal changes in the LeW and PP (Figures 2c,d and 3). In addition, the most negative seasonal mean differences fb CS2_80 -fb CS2_50 observed during these winters correspond to the highest (lowest) values of LeW (PP) which suggests a relationship between the inter-annual variability of these parameters.
The link between seasonal mean and the seasonal change rate in the LeW, PP and fb CS2_80 -fb CS2_50 is assessed by calculating correlation coefficients for each parameter ( Table 2). The correlation between the mean and the change rate for the October-to-April period is significant at the 95% confidence level for the PP, for the sea ice freeboard difference fb CS2_80 -fb CS2_50 and for the LeW derived from the CS2 data (marked with bold in Table 2). Even higher correlations are estimated between October-November Remote Sens. 2020, 12, 3094 6 of 16 mean and the change during winter seasons. Thus, we conclude that changes in surface characteristics during winter are primarily controlled by its initial conditions that are set after the preceding summer season. After intensive summer melting, a decrease in multiyear ice fraction and snow depth result in levelling sea ice surface and a decrease in volume scattering that is reflected in low LeW and high PP values. Then, the LeW rapidly increases (and PP decreases) during winter due to strong contrast between initial and subsequent surface properties as a result of surface roughness increase and snow accumulation. In the case of relatively cold summers, the surface characteristics during winter show no significant change with respect to initial conditions and changes in the waveform parameters are small.
Remote Sens. 2020, 12, x FOR PEER REVIEW 6 of 15 changes in surface characteristics during winter are primarily controlled by its initial conditions that are set after the preceding summer season. After intensive summer melting, a decrease in multiyear ice fraction and snow depth result in levelling sea ice surface and a decrease in volume scattering that is reflected in low LeW and high PP values. Then, the LeW rapidly increases (and PP decreases) during winter due to strong contrast between initial and subsequent surface properties as a result of surface roughness increase and snow accumulation. In the case of relatively cold summers, the surface characteristics during winter show no significant change with respect to initial conditions and changes in the waveform parameters are small. It should be noted that a strong relationship between the mean and the change in the surface characteristic over mooring sites is indicated despite the fact that we do not account for the sea ice drift. A spatial consistency between seasonally mean and seasonal change rate for the LeW   It should be noted that a strong relationship between the mean and the change in the surface characteristic over mooring sites is indicated despite the fact that we do not account for the sea ice drift. A spatial consistency between seasonally mean and seasonal change rate for the LeW and The seasonal mean and changes in the sea ice freeboard difference fb CS2_80 -fb CS2_50 also agree with each other, as well as with distribution of the LeW and PP waveform parameters derived from CS2 data. Thus, the consistency of spatial distributions is in agreement with the consistency of the corresponding time series for the mooring locations (Table 2). Table 2. Correlation coefficients between means (October-to-April and October-to-November) and October-to-November change rates for the LeW and PP radar altimeter waveform parameters derived from Envisat and Cryosat-2 data, for sea ice draft differences (dr Env -dr uls , dr CS2_50 -dr uls , dr CS2_80 -dr uls ) and sea ice freeboard difference fb CS2_80 -fb CS2_50 estimated using data from the mooring A and data from all mooring sites. Correlation coefficients that are significant at the 95% confidence level are denoted in bold. The seasonal mean and changes in the sea ice freeboard difference fbCS2_80-fbCS2_50 also agree with each other, as well as with distribution of the LeW and PP waveform parameters derived from CS2 data. Thus, the consistency of spatial distributions is in agreement with the consistency of the corresponding time series for the mooring locations (Table 2). Table 2. Correlation coefficients between means (October-to-April and October-to-November) and

October-April Mean
October-to-November change rates for the LeW and PP radar altimeter waveform parameters derived from Envisat and Cryosat-2 data, for sea ice draft differences (drEnv-druls, drCS2_50-druls, drCS2_80-druls) and sea ice freeboard difference fbCS2_80-fbCS2_50 estimated using data from the mooring A and data from all mooring sites. Correlation coefficients that are significant at the 95% confidence level are denoted in bold.

Impact of Sea Ice Freeboard Retrieval Uncertainly on Sea Ice Draft Estimates
Time series on Figures 2c,d and 3 show that a seasonal decrease in the differences between sea ice drafts derived from satellites and measured by ULS observed during most of the winter seasons is consistent with seasonal increase in the LeW and decrease in the PP and fbCS2_80-fbCS2_50 values. An agreement in inter-annual variability between all these characteristics can also be noted-e.g., from consistency of the low mean sea ice draft differences, especially the drCS2_80-druls, and the high (low) LeW (PP) values during 2013-2014 and 2014-2015 seasons. This suggests that sea ice freeboard retrieval uncertainties, to a large extent, can explain the seasonal and inter-annual variations of the differences between satellite-derived and ULS-measured sea ice draft. This is supported by the fact that, as well as for the LeW, PP and fbCS2_80-fbCS2_50, the correlation of seasonal changes and the mean values in the beginning of winter season are significant at the 95% confidence level for sea ice draft differences drEnv-druls and drCS2_80-druls (Table 2). Low correlation for the differences drCS2_50-druls indicates that the sea ice draft drCS2_50 is less affected by the variations in freeboard uncertainties.

Impact of Sea Ice Freeboard Retrieval Uncertainly on Sea Ice Draft Estimates
Time series on Figures 2c,d and 3 show that a seasonal decrease in the differences between sea ice drafts derived from satellites and measured by ULS observed during most of the winter seasons is consistent with seasonal increase in the LeW and decrease in the PP and fb CS2_80 -fb CS2_50 values. An agreement in inter-annual variability between all these characteristics can also be noted-e.g., from consistency of the low mean sea ice draft differences, especially the dr CS2_80 -dr uls , and the high (low) LeW (PP) values during 2013-2014 and 2014-2015 seasons. This suggests that sea ice freeboard retrieval uncertainties, to a large extent, can explain the seasonal and inter-annual variations of the differences between satellite-derived and ULS-measured sea ice draft. This is supported by the fact that, as well as for the LeW, PP and fb CS2_80 -fb CS2_50 , the correlation of seasonal changes and the mean values in the beginning of winter season are significant at the 95% confidence level for sea ice draft differences dr Env -dr uls and dr CS2_80 -dr uls ( Table 2). Low correlation for the differences dr CS2_50 -dr uls indicates that the sea ice draft dr CS2_50 is less affected by the variations in freeboard uncertainties.
To evaluate the impact of sea ice freeboard uncertainties on draft estimates, we calculated the correlation of the LeW and PP parameters derived from Envisat and CS2 data over sea ice surface as well as the sea ice freeboard difference fb CS2_80 -fb CS2_50 with the corresponding differences between Remote Sens. 2020, 12, 3094 9 of 16 satellite and ULS sea ice drafts. Table 3 presents the correlation coefficients of October to April and October to November means, as well as the October to April seasonal change rate between sea ice draft differences and the corresponding LeW, PP and fb CS2_80 -fb CS2_50 values. The highest correlation is observed for the difference dr CS2_80 -dr uls -all correlation coefficients between seasonal change rates are significant at the 95% confidence level, and notable correlations are also observed between the October and November means. For the differences dr CS2_50 -dr uls the correlation is much lower and significant only between seasonal change rates when comparing with the PP for the mooring A. For Envisat data, a sharp change in sea ice draft difference after summer 2008 is not reflected in the time series of waveform parameters, resulting in modest correlations between the dr Env -dr uls and the LeW and PP parameters, though statistically significant correlations are still observed when comparing October to April means using data from all moorings. Table 3. Correlation coefficients between differences of satellite-derived and ULS-measured sea ice drafts (dr Env -dr uls , dr CS2_50 -dr uls , dr CS2_80 -dr uls ) and the waveform parameters (LeW and PP) and sea ice freeboard differences fb CS2_80 -fb CS2_50 for October-to-April means, October-to-November means and October-to-April change rates, estimated using data from the mooring A and data from all mooring sites. Correlation coefficients that are significant are at the 95% confidence level are denoted in bold. The results presented in Tables 2 and 3 suggest that uncertainty in sea ice freeboard retrievals impacts seasonal and inter-annual discrepancies between the satellite-derived sea ice draft dr CS2_80 and the draft measured by the ULS, while the relationship for the differences between dr Env -dr uls and dr CS2_50 -dr uls is observed only in part. The found relationships can be attributed to the change in surface properties during winter that is driven by the initial conditions formed after melt season. The greater seasonal change in surface properties, as represented by the change in the difference fb CS2_80 -fb CS2_50 and waveform parameters, the larger the underestimation of seasonal sea ice draft growth. The small underestimation of sea ice draft growth for the seasons 2013-2014 and 2014-2015 corresponds to the lowest change in the LeW, PP and fb CS2_80 -fb CS2_50 values.

October-November Mean
The lower seasonal growth of the dr CS2_80 as compared to the dr CS2_50 could be attributed to the different effect of surface roughness [10] and snow accumulation [26] on elevation retrieval-the 80% retracker threshold corresponds to the lower level of scattering horizon and the estimates of ice floe elevation are less affected by the surface ridges and snow properties. This effect can explain the observed better agreement of the sea ice drafts dr uls and dr CS2_80 in the beginning of the winter seasons and with the dr CS2_50 in the end. After summer melting in the autumn months, when the sea ice in the area of the BGEP moorings is predominantly less deformed first-year ice (FYI), the better elevation estimates are provided with the retracker threshold closer to the peak of signal power, as suggested by Wingham et al. [27] for the CS2 measurements in SAR mode. In contrast, under increased roughness and snow cover in the end of winter season, the elevations are better represented by the estimates derived using the 50% retracker threshold. The dr CS2_80 is underestimated in the autumn months only in 2013 and 2014, when the LeW (PP) is high (low), indicating increased roughness during the entire winter season.
However, using the retracker threshold closer to the waveform peak is supposed to provide, on average, more accurate surface elevation estimates from the radar altimeter measurements in SAR mode over sea ice [27]. Therefore, the better agreement between the dr CS2_50 and dr uls in terms of the mean difference, rms of the difference and seasonal sea ice draft growth cannot be explained only by the uncertainties in sea ice freeboard retrieval. In particular, the comparison of satellite estimates derived using threshold-based retracking algorithms and physical model indicate that sea ice freeboard, retrieved by applying the 50% retracker threshold, is overestimated over multiyear ice (MYI), and reproduced [10] or underestimated [28] over FYI. The seasonal changes of the differences fb CS2_80 -fb CS2_50 are primarily caused by changes in the sea ice surface properties, while the application of the two different retracker thresholds for lead waveforms has a minor impact on our analysis. Lead waveforms are characterized by high PP, and echo power at the leading edge is concentrated in a few range bins at most. We, therefore, assume that choosing either a 50% or 80% threshold results in a small offset for the respective freeboard values that is constant over the winter season but does not change the divergence between the two freeboard estimates that is primarily controlled by the increasing LeW over ice surfaces.
To evaluate the effect of the sea ice type on the draft estimates derived using different retracker thresholds, we used the data on the MYI fraction, provided by the Ocean and Sea Ice Satellite Application Facility, and for the SIT retrieval in the course of satellite altimeter data processing. Figure 6 shows, as expected, that time series of the MYI fraction over BGEP mooring sites increases during most of the winter seasons, providing information about the dependence of the sea ice draft estimates on the sea ice type. In contrast to the analysis of the sea ice freeboard reported in [10,28], the sea ice draft difference dr CS2_50 -dr uls indicates that the satellite-derived draft is overestimated in the beginning of winter season and becomes closer to the ULS-measured draft with the increase in MYI fraction.
represented by the estimates derived using the 50% retracker threshold. The drCS2_80 is underestimated in the autumn months only in 2013 and 2014, when the LeW (PP) is high (low), indicating increased roughness during the entire winter season.
However, using the retracker threshold closer to the waveform peak is supposed to provide, on average, more accurate surface elevation estimates from the radar altimeter measurements in SAR mode over sea ice [27]. Therefore, the better agreement between the drCS2_50 and druls in terms of the mean difference, rms of the difference and seasonal sea ice draft growth cannot be explained only by the uncertainties in sea ice freeboard retrieval. In particular, the comparison of satellite estimates derived using threshold-based retracking algorithms and physical model indicate that sea ice freeboard, retrieved by applying the 50% retracker threshold, is overestimated over multiyear ice (MYI), and reproduced [10] or underestimated [28] over FYI. The seasonal changes of the differences fbCS2_80-fbCS2_50 are primarily caused by changes in the sea ice surface properties, while the application of the two different retracker thresholds for lead waveforms has a minor impact on our analysis. Lead waveforms are characterized by high PP, and echo power at the leading edge is concentrated in a few range bins at most. We, therefore, assume that choosing either a 50% or 80% threshold results in a small offset for the respective freeboard values that is constant over the winter season but does not change the divergence between the two freeboard estimates that is primarily controlled by the increasing LeW over ice surfaces.
To evaluate the effect of the sea ice type on the draft estimates derived using different retracker thresholds, we used the data on the MYI fraction, provided by the Ocean and Sea Ice Satellite Application Facility, and for the SIT retrieval in the course of satellite altimeter data processing. Figure 6 shows, as expected, that time series of the MYI fraction over BGEP mooring sites increases during most of the winter seasons, providing information about the dependence of the sea ice draft estimates on the sea ice type. In contrast to the analysis of the sea ice freeboard reported in [10,28], the sea ice draft difference drCS2_50-druls indicates that the satellite-derived draft is overestimated in the beginning of winter season and becomes closer to the ULS-measured draft with the increase in MYI fraction.
Another possible reason for the observed underestimation of seasonal growth of the SIT derived from satellite radar altimetry can be associated with the fact that the underestimation of seasonal growth of the SIT derived from satellite radar altimetry is also associated with incorrect accounting for the lower speed of radar wave propagation in snow that is applied in the course of freeboard retrieval [29]. However, this effect is too small to explain the observed differences between satellite and ULS sea ice draft.  Another possible reason for the observed underestimation of seasonal growth of the SIT derived from satellite radar altimetry can be associated with the fact that the underestimation of seasonal growth of the SIT derived from satellite radar altimetry is also associated with incorrect accounting for the lower speed of radar wave propagation in snow that is applied in the course of freeboard retrieval [29]. However, this effect is too small to explain the observed differences between satellite and ULS sea ice draft.

Impact of Sea Ice Freeboard to Thickness Conversion on Sea Ice Draft Estimates
Based on comparison of the dr CS2_50 and dr CS2_80 time series with the ULS data, we argue that the uncertainties related to the conversion of freeboard to thickness should significantly impact the satellite-derived draft estimates. The thickness (h ice ) is calculated from the hydrostatic equilibrium equation: where h f r is sea ice freeboard, h s is snow depth, and ρ w , ρ ice and ρ s are the densities of water, ice and snow. The largest uncertainties in the SIT estimates, along with those related to the sea ice freeboard retrieval, are associated with poor knowledge about snow depth and sea ice density [30,31]. In principle, the underestimation of the sea ice draft seasonal growth can be attributed to an underestimation of the seasonal snow depth growth by Warren climatology [11], used in the hydrostatic equilibrium equation to account for the snow loading and potentially changes in bulk sea-ice density.
The largest underestimation of the sea ice draft seasonal growth is observed for the dr CS2_80 that correlates with seasonal changes of the LeW and PP and fb CS2_80 -fb CS2_50 parameters ( Table 3). Assuming that the underestimation of the sea ice draft seasonal growth is caused by the underestimation of snow depth growth in climatological data, the seasonal changes in the waveform parameters and the difference in fb CS2_80 -fb CS2_50 should also, apparently, result from snow accumulation during winter. However, in this case, reasonably reproduced seasonal sea ice draft growth during winters 2013-2014 and 2014-2015 could not result from low snow accumulation indicated by the observed small seasonal change in the LeW, PP and fb CS2_80 -fb CS2_50 . Moreover, it is unlikely that Warren climatology underestimates the snow depth growth for most of the winter seasons over the considered time-periods of satellite and ULS observations.
The uncertainties in the sea ice density result from using fixed values representing typical densities of the FYI and MYI [12]. Belter et al. [15] showed that the strongest underestimation of the sea ice draft derived from the CS2 data corresponds to an increase in the sea ice deformation. Assuming that the ULS measures the distance to the underside of the blocky sea ice structures, one needs to consider the significant sea-water filled pore spaces in the density as well. In contrast, consolidated multi-year ice ridges may have a different bulk density closer to the level ice densities than their young unconsolidated counterparts. Therefore, since sea ice deformation increases during the winter seasons, as indicated from changes in waveform parameters, this may lead to a seasonal increase in the underestimation of sea ice density and, hence, SIT estimates. This effect could also explain inter-annual changes in the difference between satellite-derived and ULS-measured draft. In particular, the sharp change in sea ice draft difference after summer 2008 could be due to the systematic underestimation of sea ice density of heavily deformed ice during preceding years, though time series of the waveform parameters derived from Envisat measurements do not reflect significant change in surface properties. It should be noted that, in contrast to snow depth, which also affects the sea-ice freeboard estimate, sea-ice density is a parameter solely relevant in the freeboard to thickness conversion. Therefore, better agreement between the dr uls and dr CS2_50 can be attributed to the systematic underestimating of sea ice density and its change over the winter season combined with the overestimating of the freeboard fr CS2_50 with the two error contribution partly compensating each other.

Effect of Summer Conditions and Sea Ice Type on Sea Ice Draft Estimates
To investigate the initial conditions that precede winter seasons we examined monthly-averaged sea surface temperature (SST) in the area of the BGEP moorings obtained from the ERA-5 reanalysis data provided by the European Centre for Medium-Range Weather Forecasts for the 0.25 • × 0.25 • grids. For comparison with the difference between satellite retrievals and ULS measurements, we used the SST data averaged within the radius of 200 km around mooring locations.
When the ocean is covered by sea ice, the SST in ERA-5 reanalysis is set to 271.46K and, therefore, does not provide any valuable information for winter months. In summer, at least part of the area surrounding BGEP moorings have open water, and SST data can be used as a measure of melt intensity ( Figure 6). To characterize the sea ice surface during winter season, we also used time series of MYI fraction that, in particular, strongly depends on melt intensity in the preceding summer.  Figure 6 also shows that the sea ice type during winter seasons following summers 2008 and 2013 is predominantly FYI and MYI, respectively. Table 4 presents correlation coefficients of the differences between satellite and ULS sea ice drafts with July and August SST and with mean and change in MYI fraction during winter season. In addition to correlations for the three considered satellite-derived sea ice draft estimates, Table 4 presents correlations for the combined time series derived from inter-mission consistent sea ice drafts, the dr Env and dr CS2_50 , by averaging values during the overlap period. The sea ice draft differences after summer 2008 remained relatively high for the rest of the dr Env -dr uls and for the whole dr CS2_50 -dr uls time series, which may indicate some persistent change in sea ice conditions influencing the results of the sea ice thickness retrieval procedure. Since it may prevent evaluation of the relationships for the entire time series the correlations are calculated, in addition, for time series formed from the differences between adjacent points of the original time series. Similar to the comparison with the LeW, PP, and fb CS2_80 -fb CS2_50 parameters, the relationship is most consistent for the sea ice draft difference dr CS2_80 -dr uls , and the largest correlation coefficients are observed for October-November mean differences. Correlation for the dr Env -dr uls is similar for October-November and October-April mean values, but it is statistically significant only for time series formed from the differences between adjacent points. The notable effect of using this type of time series is seen also for the combined sea ice draft time series that provides a significant correlation at 95% confidence level for the mean draft differences. The lowest relationship is estimated for the difference dr CS2_50 -dr uls with notable correlation observed only between their seasonally mean values and summer SST when combining data from all moorings.
The correlation of seasonal changes is lower for all differences between satellite and ULS sea ice drafts, although limited spatial consistency between mean values and seasonal change rates still can be  The correlation of seasonal changes is lower for all differences between satellite and ULS sea ice drafts, although limited spatial consistency between mean values and seasonal change rates still can be noted for the MYI

Conclusions
In this study, we used the BGEP ULS measurements in the Beaufort Sea for the period 2003-2018 to validate sea ice draft time series derived from the Envisat and CS2 satellite radar altimeter data by the algorithm used to produce the ESA CCI-2 SIT climate data record. The effect of sea ice  Table 4. Correlation coefficients between differences of satellite-derived and ULS-measured sea ice drafts (dr Env -dr uls , dr CS2_50 -dr uls and dr CS2_80 -dr uls ) as well as combined time series derived from the differences between dr Env -dr uls and dr CS2_50 -dr uls , and the MYI fraction and SST in July and August for October-to-April means, October-to-November means and October-to-April change rates, estimated using data from the mooring A and data from all mooring sites. Correlation coefficients presented in columns 1 and 2 are derived from original time series and from time series formed using the differences between adjacent points of the original time series. Correlation coefficients significant at the 95% confidence level are denoted in bold.

Conclusions
In this study, we used the BGEP ULS measurements in the Beaufort Sea for the period 2003-2018 to validate sea ice draft time series derived from the Envisat and CS2 satellite radar altimeter data by the algorithm used to produce the ESA CCI-2 SIT climate data record. The effect of sea ice freeboard uncertainties due to changes in surface properties on the sea ice draft satellite estimates was analysed using time series of the waveform parameters and the difference between sea ice freeboards retrieved from CS2 data by applying 50%, as used for the CCI-2 SIT product, and 80% retracker thresholds. Time series indicate that all satellite retrievals underestimated seasonal sea ice draft growth during most of winter seasons over the considered time period. The seasonal mean and seasonal changes in the difference between satellite-derived and ULS-measured sea ice draft are found to be linked with the sea ice surface properties preconditioned by the melt intensity during the preceding summer. The best agreement with the ULS data in terms of bias, rms of the differences and underestimation of sea ice draft growth during winter is indicated for the results derived using the 50% retracker threshold of the existing CCI climate data record. Since more accurate elevation estimates from the CS2 measurements in SAR mode are supposed to be provided by the 80% retracker threshold, this result can be attributed to the combined effect of overestimating sea ice freeboard and underestimating sea ice density used for the conversion of freeboard to thickness on SIT estimates.
Sea ice retrievals from the satellite radar altimetry using different algorithms, based on the threshold-based retracker approach [5,19,23] and physical model [4,10,28], provide a range of sea ice freeboard and thickness estimates-e.g., [20,32]. The impact of surface properties on satellite-derived sea ice retrievals and the analysis of its origin presented in this study may be used for the evaluation and future development of retrieval algorithms.