Climate change is among the most pressing issues of our time, with many severe economic, environmental, and geopolitical consequences. Recently, the application of time series analytical methods to this topic – and, more broadly, a “climate econometrics” – has emerged as a vibrant research literature, as highlighted, for example, in Hillebrandetal2020 and the references therein. One important issue that these methods can address is the diminishing area or extent of Arctic sea ice. The loss of Arctic sea ice is a vital focus point of climate study because it is both an ongoing conspicuous effect of climate change and a cause of additional climate change via feedback loops, including the fact that reduced Arctic sea ice boosts solar energy absorption (DRice).
There are, however, several alternative measures of Arctic sea ice extent based on different processing methodologies of the underlying satellite-based microwave measurement data, and the choice among these measures is not clear-cut (BunzelEtAl16). In this paper, we study four such sea ice extent () measures, which we denote as Sea Ice Index (), Goddard Bootstrap (), JAXA (), and Bremen (). The top panel of Figure 1 provides time series plots of these four measures of Arctic for the satellite measurement era, which started in 1978. The four measures appear almost identical, because their scale is dominated by large seasonal swings. However, the effects of seasonality can be removed by plotting each month separately for the four series, as done in the lower twelve panels of Figure 1. Of course, the Arctic measures all trend down in every month, with steeper trends for the summer months (e.g., August, September, October). (Note the different axis scales for different months.) There are also systematic differences across indicators. , for example, tends to be high, and tends to be low, while and are intermediate. But the deviations between various pairs of measures are not rigid; that is, they are not simply parallel translations of each other. Instead there are sizable time-varying differences among the various measures.
All of this suggests treating the various measures as noisy “indicators” of a latent true sea ice extent, which in turn suggests the possibility of blending them into a single combined indicator with less measurement error. Indeed some prominent studies have used simple equally-weighted averages of competing indicators, with precisely that goal. For example, a recent report on the state of the cryosphere (IPCC19) uses a simple average of three indicators.111See the notes for Figure 3.3, page 3-13. Simple averages, however, are often sub-optimal. Optimality generally requires use of weighted averages giving, for example, less weight to noisier indicators. Motivated by these considerations, in this paper, we propose and explore a dynamic factor state-space model that combines the various published measures into an optimal measure of sea ice extent, which we extract using the Kalman smoother.
We proceed as follows. In section 2, we describe the four leading Arctic sea ice extent indicators that we study, and the satellites, sensors, and algorithms used to produce them. In section 3, we propose a basic dynamic-factor state-space model for sea ice extent and use it to obtain optimal extractions of latent extent. We conclude in section 4.
2 Four Arctic Sea Ice Extent Indicators
Sea ice extent () indicators are constructed from satellite measurements of the earth’s surface using passive microwave sensing, which is unaffected by cloud cover or a lack of sunlight. Several steps are necessary to convert raw reflectivity observations into final measurements. First, for a polar region divided into a grid of individual cells, various sensors record a brightness reading or “brightness temperature” for each cell. An algorithm then transforms these brightness readings into fractional surface coverage estimates – sea ice concentration () values – for each grid cell. Finally, is calculated by summing the area of all cells with at least 15 percent ice surface coverage.222PARKINSONETAL2008 discusses reasons for using a 15 percent cutoff. This up-rounding in is effectively a bias correction, as determining the edge between ice and water can be especially difficult in the summer, when, for example, melting pools on summer ice surfaces can be mistaken for ice-free open water (MeierAndStewart2019).
Different algorithms for processing the raw measurements importantly shape the final estimates. In addition, the seiries are based on raw data that can be obtained from somewhat different satellites and sensors (ComisoEtAl17, NSIDC_Boot_Explained). In this section, we review some aspects of the satellites, sensors, and algorithms that underlie the measures.
2.1 Satellites and Sensors
Table 1 summarizes the operative dates of the various satellites and sensors relevant for Arctic ice measurement. The first multi-frequency sensor equipped on a satellite was the Scanning Multichannel Microwave Radiometer (SMMR) launched in 1978 (NSIDC_Caval). Starting in 1987, later sensors – the Special Sensor Microwave Imager (SSM/I) and the Special Sensor Microwave Imager/Sounder (SSMIS) – offered higher resolution images.333For detailed discussion of sensor characteristics see https://nsidc.org/ancillary-pages/smmr-ssmi-ssmis-sensors. In 2002 and 2012, respectively, the Advanced Microwave Scanner Radiometer for EOS (AMSR-E) and Advanced Microwave Scanner Radiometer 2 (AMSR2) sensors were launched and provided further improvements in resolution (ComisoEtAl17).444Early in the sample, operational problems prevented data delivery for several days during 1986 and between December 1987 and January 1988 (NSIDC_Boot). For more recent technical difficulties, see https://www.nrl.navy.mil/WindSat/Description.php. Given the inclinations of the satellite orbits and the spherical shape of the earth, all of the satellites share an inability to observe the Arctic pole hole – a circular region at the very top of the world. The size of the pole hole varies across sensors, but historically, there is full confidence that it fulfills the 15 percent requirement (MeierAndStewart2019).
Table 1 also describe the underlying source data for our four indicators. These measures use algorithms to transform the raw satellite brightness data into and values. We now turn to a more detailed discussion of these algorithms to illuminate the differences across indicators.
|Satellite / Sensor||NASA Team||Goddard Bootstrap||JAXA Bootstrap||ASI|
Notes: NASA Team and Goddard Bootstrap dates from NSIDC_SII. JAXA Bootstrap dates from https://kuroshio.eorc.jaxa.jp/JASMES/climate/index.html. ASI dates from https://seaice.uni-bremen.de/sea-ice-concentration/time-series/.
2.2 From Raw Measurement to : Algorithmic Transformations
Once brightness data have been recorded by satellite sensors, an algorithm converts the measurements into estimates of . Here, we discuss the algorithms and other details of the various indicators.
2.2.1 Sea Ice Index
Updated on a daily basis and distributed by the National Snow & Ice Data Center (NSIDC), the Sea Ice Index (SII or ) combines two separate Sea Ice indicators: (1) the Sea Ice Concentrations from Nimbus-7 SMMR and DMSP SSM/I-SSMIS Passive Microwave Data (NASA Team) (NSIDC_Caval) – produced at the Goddard Space Flight Center – and (2) the Near-Real-Time DMSP SSMIS Daily Polar Gridded Sea Ice Concentrations (NRTSI) (NSIDC_NearReal) – produced by the NSIDC itself.555See https://doi.org/10.7265/N5K072F8. A time-lag of about one year between the estimates by NASA Team and its publication in the NSIDC database requires the NRTSI to complement the SII.
The NRTSI follows the NASA Team algorithm as closely as possible, but inconsistencies between the two series cannot be ruled out entirely (NSIDC_SII). In particular, the two sub-indicators use brightness temperatures from different providers.666NSIDC_NearReal takes the data from the National Oceanic and Atmospheric Administration Comprehensive Large Array-data Stewardship System (NOAA CLASS); (NSIDC_Caval) uses data processed at the NASA Goddard Space Flight Center. These raw readings can be distorted by weather effects, making open water look like sea ice cover. Therefore, post-calculation quality checks apply land and ocean masks, to remove errorneous and implausible ice covers. However, NASA Team and NRTSI do not apply the exact same filters (NSIDC_SII). The former algorithm additionally screens the data manually for falsely detected ice formation (NSIDC_Caval), which can enhance accuracy but also reduce transparency of the final measurements.
As Table 1 shows, the SII obtains the raw data from different generations of satellites and sensors. To make the data comparable, a linear least-squares model on the brightness temperatures, as reported by the two distinct sensors for an overlapping period of operation, is intended to adjust the reference points of 100 percent sea ice and 100 percent open water. These tie points then remain fixed over the lifetime of the new system (NSIDC_Caval_tie).
2.2.2 Goddard Bootstrap
Another sea ice indicator, distributed by the NSIDC, relies on estimates of the Goddard Bootstrap algorithm ().777For detailed algorithm description see ComisoEtAl17b. For the data, see https://doi.org/10.5067/7Q8HCCWS4I0R. Despite the NASA Team and the Goddard Bootstrap algorithms having both been developed at the NASA Goddard Space Flight Center, there are some differences between the two approaches. These arise mostly from the calibration of tie points: Other than the NASA Team, which adjusts these reference points for 100 percent open water and 100 percent ice only when a new satellite or sensor becomes operational, the Goddard Bootstrap algorithm adjusts these reference points on a daily basis to account for varying weather conditions (ComisoEtAl17). Differing weather filters, and sensitivities to varying physical temperature also lead to differences in the final measurements(ComisoEtAl97). In contrast to the NASA Team, the strength of the Goddard Bootstrap algorithm is the identification of melting sea ice. Therefore, the Goddard Bootstrap algorithm provides more accurate estimates of the edge of the ice cover (goldsteinetal2018).
Although differences between and are generally assessed to be small, they cannot necessarily be neglected (goldsteinetal2018). The differences between the NASA Team and Goddard Bootstrap algorithms occur especially during the melting season, when the former generally reports larger deviations from ship or radar observations. However, the relative accuracy of the two algorithms is not clear-cut, as the Goddard Bootstrap algorithm is highly sensitive to physical temperature and underestimates during winter times in the higher latitudes of the Arctic region.
2.2.3 Japan Aerospace Exploration Agency
As listed in Table 1, both and rely on the same set of instruments, which have been criticized for their low spatial resolution (goldsteinetal2018). The Japan Aerospace Exploration Agency (JAXA) sea ice measure () uses an adapted version of the Goddard Bootstrap algorithm to derive measures from satellite readings with higher spatial resolution (ComisoEtAl17b) .888For description and data see https://kuroshio.eorc.jaxa.jp/JASMES/climate/index.html. However, readings from these high resolution satellites are only available since 2000, so their data must merged with observations from older sensors to extend the data coverage to 1978 (ComisoEtAl17). also distinguishes itself from and by using a moving average of five days of observation to compensate for potentially missing data.
2.2.4 University of Bremen
Using observations delivered by the high-resolution AMSR-E sensor, a group of researchers at the University of Bremen developed the ARTIST Sea Ice (ASI) algorithm (SpreenEtAl08) to estimate daily .999Monthly data from https://seaice.uni-bremen.de/data/amsr2/today/extent_n_19720101-20181231_amsr2.txt. The time-series () uses different algorithms for different sensors. Until the launch of the AMSR-E sensor in 2003, used the NASA Team algorithm to transform brightness readings into values. From then on, the NASA Team algorithm is replaced by the ASI.101010See https://seaice.uni-bremen.de/sea-ice-concentration-amsr-eamsr2/time-series/.
3 Optimal Extraction of Latent Extent
The four sea ice measures discussed above differ in terms of the raw data sources and the algorithms used to process those data. They can be viewed as distinct indicators of an unobserved or latent “true” sea ice extent, , and blending such a set of noisy indicators into a composite can produce a single series with less measurement error. Here we formalize this intuition in a state-space dynamic-factor model, from which we extract an optimal composite estimate of from the four component indicators.
3.1 A Dynamic Factor Model
We work in a state-space environment, modeling each of the four indicators (, , , and ) as driven by latent true sea ice extent, , with an additive measurement error.111111The approach parallels ADNSS2016, who extract latent U.S. GDP from noisy expenditure-side and income-side estimates.
The measurement equation is
The transition equation is
where is orthogonal to at all leads and lags. Various modeling approaches are distinguished by their treatment of and . We follow DRice and allow for 12 monthly deterministic seasonal effects, each of which is endowed with (possible) deterministic quadratic trend.121212Alternatively, one could entertain stochastic trend and seasonality, but we leave that to future research. A simple approach would be separate month-by month modeling so that there is no seasonality, whether with one unit root, as in (for month ) , or two unit roots, as in , where . This results in a blended deterministic “trend/seasonal” given by
where indicates month and indicates time. Hence the full transition equation is
The model is already in state-space form, and one pass of the Kalman filter, initialized with the unconditional state mean and covariance matrix, provides the 1-step prediction errors necessary to construct the Gaussian likelihood. We maximize the likelihood using the EM algorithm and calculate standard errors from the analytic Hessian matrix. Following estimation, we use the Kalman smoother to obtain the minimum variance unbiased extraction offrom the estimated model. The smoother averages across indicators, but it desirably produces optimally weighted averages rather than simple averages. The smoother also averages over time, using data both before and after time to estimate , which is also necessary for minimum-variance extraction, due to the serial correlation in . For details see Harvey1989.
One or more restrictions are necessary for identification. The standard approach is to normalize a factor loading, which amounts to an unbiasedness assumption. Normalizing , for example, amounts to an assumption that is unbiased for . Whether there truly exists such an unbiased indicator (and if so, which) is of course an open question – one can never know for sure. and are the most widely used indicators (MeierEtAl14, PengEtAl13), so it is natural to consider normalizing on or . We explore both.
3.2 Estimated Measurement Equation
The estimated measurement equation (1), normalized with , is
where standard errors appear beneath each estimated loading. All indicators are estimated to load heavily on , with all very close to 1. loads least heavily (), in accord with its generally lower values in Figure 1, and loads most heavily (), in accord with its generally higher values. loads with an estimated coefficient insignificantly different from 1, and, of course, the loading is 1 by construction. Hence, the estimation results accord with Figure 1, with and more in the center of the range and and more extreme.
Alternatively, the estimated measurement equation normalized with is
The estimated loadings in equations (7) and (8), corresponding to and respectively, are the same up to the normalization. All estimated loadings are now less than one, because we now normalize on Goddard, which is always the highest.
Now consider the associated measurement error covariance matrix (3). The estimate for the normalization is
with implied estimated correlation matrix
Note that is much higher than any of , , and , potentially due to different indicators using different methods to determine tie points, i.e., reference points of brightness for 100% sea ice and 100% open water. The choice is crucial for accurate measurement of within grid cells. Tie points, moreover, need not be constant, as brightness readings are sensitive to weather effects and atmospheric forcings (Ivan15). Dynamic tie-point calibration is potentially desirable because it can decrease the bias of measurements (ComisoEtAl17). The latest version of the Goddard Bootstrap algorithm, in particular, calibrates tie points daily. One would expect, however, that the bias reduction from dynamic tie-point calibration may come at the cost of potential discontinuities that increase measurement error variance. Our results confirm that conjecture. Our estimate of is about 4.5 times that of (which uses constant tie points).
Alternatively, the estimated measurement error covariance matrix for the normalization is
with implied estimated correlation matrix
3.3 Estimated Transition Equation
Notes: The are and the are .
Now let us move to the transition equation (6). Using the normalization we obtain and trend/seasonal parameter estimates (, , and ) as reported in the columns of Table 2. Trends for all months are highly significant and downward sloping. The trends for summer months (August, September, October, November) display a notable negative, and generally statistically significant, curvature, whereas for non-summer months, the quadratic trend terms are generally small and statistically insignificant. For the normalization we get and trend/seasonal parameter estimates (, , and ) as reported in the columns of Table 2. The and results are very similar.
3.4 Extracted Latent Sea Ice Extent
In Figures 2 () and 3 () we show we show optimal latent sea ice extent extractions () in black, together with the four raw indicators in color, by month. First consider Figure 2. Of course is centered on due to the normalization. Moreover, is always very close – almost identical – to (and close to , because tends to be very close to ).
Now consider Figure 3 (). Due to the different normalization, is centered not on but rather on , so is shifted upward relative to . The location of relative to , moreover, clearly varies by month. In winter months, it tends to be greater than , whereas in summer months it tends to be less than . Note in particular that the variation of around is noticeably greater than the variation of around . Clearly is influenced more by movements in other indicators (, , and ) than is .
The fact that and are different, both in terms of level and variation around the level, limits their usefulness for research focusing on level, as the level depends entirely on identifying assumptions. However, and crucially, in an important sense and are highly similar: The model and identification scheme makes
and identical up to a linear transformation.
identical up to a linear transformation.This is clear in Figure 4, which plots the two competing extracted factors, month-by-month. Regressions of on yield highly-significant intercepts not far from zero, highly-significant slopes near 1.0, and values above 0.999, for each month. The amount of the intercept shift varies by month, and is larger in winter than summer. In addition, the difference between the two extracted extents trends downward in all months, as shown in Figure 5.
Because and are identical up to a linear transformation, it makes no difference which we use for research focused on linear relationships between and other aspects of climate (e.g., various radiative forcings). The obvious choice, then, is , which is just itself, dispensing with the need to estimate the factor model.
4 Summary and Conclusion
We propose a dynamic factor model for several Arctic sea ice extent indicators (Goddard, Sea Ice Index, JAXA, Bremen). We estimate the model and use it in conjunction with the Kalman smoother to produce a statistically-optimal combination of the individual indicators, effectively ”averaging out” the individual measurement errors. We explore two identification strategies corresponding to two different factor loading normalizations. The corresponding two extracted combined measures (latent factors) are identical up to a linear transformation, so either one can be used to explore relationships of sea ice extent and other variables. Interestingly, however, the extracted factor for one of the normalizations puts all weight on the Sea Ice Index. Hence the Sea Ice Index alone is a statistically optimal “combination” and one can simply use it alone with no loss, dispensing with the need to estimate the factor model. That is, there is no gain from combining the Sea Ice Index with other indicators, confirming and enhancing confidence in the Sea Ice Index and the NASA Team algorithm on which it is based, and similarly lending credibility – in a competition against very sophisticated opponents – to the NSIDC’s claim that the Sea Ice Index is the “final authoritative SMMR, SSM/I, and SSMIS passive microwave sea ice concentration record” (NSIDC_SII).