Glacial Cover Assessment from 1993 to 2021 Using Normalized Difference Snow Index with Landsat Imagery in the Colorado Front Range, USA ()
1. Introduction
Throughout the world, glaciers are retreating rapidly, and the mass losses have also increased in recent decades (Gardner et al., 2013) . Some glaciers are in danger of disappearing altogether in the coming decades (Viola, 2019) . Glacier ice covers about 10% of Earth’s surface and 68.7% of freshwater storage (Gleick & White, 1993) . Changes in the volume and area of the glaciers are important indicators of global change (Barry, 1999) . Some warm-based cirque glaciers may also be a substantial contributor to the water supply at local, regional, and global scales (Huss, 2011; Jansson et al., 2003; Bahr & Meier, 2000) . Moreover, glacier change affects global sea level rise and land water resources. For these reasons, it is extremely important to document the changes in glacier cover (Fountain et al., 1997) .
More than 60 percent of Colorado’s annual precipitation falls as snow, and snowmelt runoff produces about 80 percent of streams (Fassnacht et al., 2018; Serreze et al., 1999) . Since the 1990s, the shrinkage of cirque glaciers at Rocky Mountain National Park has been occurring rapidly and most of the glaciers have exhibited a significant recession since 1999 (Hoffman et al., 2007) . It was found that due to the increasing temperature, the cirque glaciers had a more-or-less shrinkage during the last century (Dyurgerov & Meier, 2000) . Satellite images of glacial areas provide one method of monitoring glacial activity in the region. By satellite remote sensing, direct information about the physically inaccessible surface resources can be obtained (Fountain et al., 1997; Papa et al., 2013) . Landsat observations are readily available and could be appropriate for cirque glaciers, although the coarser resolution could result in less accurate measurements. However, data from Landsat are widely used to improve and validate snow cover rate inversion (Appel, 2018) . Moreover, especially in the past 30 years, the application of satellite remote sensing has made generational progress in monitoring surface resources (Alsdorf & Lettenmaier, 2003; Smith, 1997) .
The advent of remote sensing sciences has opened up new opportunities for glacier studies. With the advancement of this field, satellite images now may have high resolution enough to show small glaciers existing in difficult topographic terrain. Due to the complexity of the local terrain and the relatively small size of cirque glaciers, these are generally considered unsuitable for analysis using coarser-resolution satellites (Kuhn, 1995) . Because of the long experience with usage, and scientific advancement, the resolution of satellite images is improving (Fountain et al., 1997) . In this study, first, the feasibility of using Landsat to track cirque glacier areas using late summer imagery has been analyzed. In addition, the Normalized Difference Snow Index (NDSI) has also been calculated to distinguish snow from other land cover types (Riggs et al., 2006) in the late summer when the greatest area of snow-free ice is exposed (Fountain et al., 1997) and the trend from 1993 to 2021 has also been analyzed.
2. Materials and Methods
The Front Range lies on the southeastern most region of the Rocky Mountains in the west-central United States, ranging from the southeastern Wyoming state to the south-central of Colorado State with 65 to 80 kilometers width. It has a large elevation gradient and a large amount of seasonal snowpack (Kampf & Lefsky, 2016) . There are 14 named glaciers present in the Front Range (Figure 1), many located in Rocky Mountain National Park: Arapaho Glacier, Isabelle Glacier, Dove Glacier, Mills Glacier, Sprague Glacier, Taylor Glacier, Andrews Glacier, Fair Glacier, Moomaw Glacier, Peck Glacier, Rowe Glacier, Saint Vrain Glacier, Tyndall Glacier and Saint Mary Glacier. These glaciers occur in east- to north-facing cirques at higher altitude near the ridge tops. These topographic settings favors snow accumulation from avalanches and windblown snow in the winter and are well-shaded from the summer sun. Surface air temperature during the summer and snow accumulation in the winter mostly determines glacier surface mass balance (Larsen et al., 2007) .
Of the 14 named glaciers present in the Front Range, 13 glaciers have been considered for this study. To study the detailed analysis of glacier area, satellite imageries were selected from the year 1993 to 2021. During the selection of images, it was ensured that 1) cloud cover and haze were minimal or absent, 2) fresh snow was limited, and 3) the images fell within a month of all others. In this study, the NDSI used to calculate snow cover area is computed by the following equation:
Figure 1. Satellite image of the Colorado Front Range Region, showing area under study (map source-Google Earth Landsat/Copernicus).
The following instruments from Landsat 5, Landsat 7, and Landsat 8 were used to create the multispectral satellite imagery and to calculate the NDSI. Details on the Landsat satellites instruments are included in the table below:
Thematic Mapper (TM) on Landsat 5.30 m spatial resolution of visible and near-infrared bands in operation since 1984 and decommissioned in 2013. Landsat 5 delivered Earth-imaging data for nearly 29 years.
Enhanced Thematic Mapper Plus (ETM+) on Landsat 7.30 m spatial resolution of visible, near-infrared, and short-wave infrared (SWIR) bands in operation since 1999.
Operational Land Imager (OLI) on Landsat 8.30 m spatial resolution of visible, near-infrared, and shortwave infrared (SWIR) in operation since 2013.
In order to classify the image as glacier versus not glacier, we used a threshold NDSI of 0.4 and this is the global reference NDSI threshold. For the accuracy of the observations, the glacier area calculated from Landsat multispectral satellite imagery was compared with a 41cm high-resolution Worldview satellite powered by MAXAR on the same date for all of the 13 glaciers in this study. I digitized the visible glacier area in the high-resolution Worldview images and calculated the area. With the Landsat images, I also calculated the area of 13 glaciers in the late summer when the greatest area of snow-free ice is exposed from 1993 to 2021. The glacier area of 13 glaciers was summed and tested for a significant trend over time using a Mann-Kendall test. Pearson correlation between the total area of 13 glaciers and the precipitation and annual maximum temperatures was also tested. The precipitation data and the annual maximum temperature used is from the Grand Lake 6 SSW Weather Station (40˚10'59.87"N, 105˚52'2.06"W), operated by the USGS. Grand Lake 6 SSW Weather Station is located 17.01 kilometers west of Saint Vrain Glacier, 25.64 kilometers northwest of Arapaho Glacier, and 19.34 kilometers southwest of Andrews Glacier.
Mann-Kendall trend test is based on the correlation between the ranks of a time series and their time order (Hamed, 2008) . The distribution-free Mann-Kendall test is widely used for the assessment of the significant trends in many hydrologic and climatic time series (Hamed, 2009) . This trend test has been applied in this study to understand the trend in the areal change of the glaciers.
3. Results and Discussion
To understand the feasibility of using Landsat images for this study, a 41 cm high-resolution image powered by MAXAR in Google Earth Pro has been obtained and compared with the image on the same day from Landsat (Earth, 2021). The results have shown that Landsat can track the cirque glaciers relatively well. The glacier area calculated from Landsat multispectral satellite imagery was compared with the 41 cm high-resolution Worldview satellite powered by MAXAR on the same date. It provided a weighted arithmetic mean percent error of 13 cirque glaciers in the Front Range region of 9.98% (Table 1, Figure 2). In general, the smaller glaciers, such as Mills glacier, Peck glacier, and Dove glacier showed a relatively higher percent error. Larger glaciers, such as the Saint Vrain glacier and Fair glacier showed a relatively smaller percent error (Table 1, Figure 2). Thus, it appears possible to calculate the size of cirque
Table 1. The percent error of 13 glaciers in the Colorado Front Range region on July 31, 2016.
Figure 2. Graphical representation of the calculated percent error of the area of 13 glaciers in the Colorado Front Range region on July 31, 2016. Black bar is the area calculated from 30 m-resolution Landsat images on July 31, 2016. Grey bar is the area calculated from 41 cm-resolution Worldview satellite images on same date.
glacier by using Landsat satellite imagery.
The Mann-Kendall trend test showed that the summed total area of 13 glaciers in the Colorado Front Range region is significantly decreasing (p-value < 0.05) (Table 2). Kendall rank correlation coefficient (tau) for the analysis was found 0.266 (p-value < 0.05), which indicates a high correlation.
Table 2. A significant decreasing trend of the total area of the glaciers in Front Range, USA.
The changes in glacier area in the Colorado Front Range region broadly follow a linear regression model (Hoffman et al., 2007) . The linear regression of the total area of 13 glaciers in the Colorado Front Range region is shown by the following equation having R2 value of 0.1881.
The equation shows the total area of 13 glaciers has decreased by 63.64% from 1993 to 2021 (Figure 2). The negative slope also suggests a negative correlation between area of the glaciers and years. Thus, with increasing time the glacial spatial cover was found to be decreasing, which is shown in the Figure 2. The erratic distribution of glacial cover during different years is due to change in seasonal patterns such as winter prevailing winds have a strong westerly component that redistributes snow from broad, unglaciated peneplain surfaces above treeline into the east-facing cirques (Outcalt, 1965) .
A significant correlation among the trend of the glacier area to the decreasing precipitation and increasing maximum temperatures has been detected (Figure 3). In a previous study also, the precipitation data and annual maximum temperature data of Grand Lake 6 SSW Weather Station were highly correlated with the area of the glaciers (Hoffman et al., 2007) . The total area of the glacial cover showed a negative correlation with the annual maximum temperature. In 1993 when the annual maximum temperature was below 50 (F), the total spatial glacial cover was above 2 km2. Total glacial area showed some peaks in the year of 1996, 1998, 2008, and 2011, which were relatively colder years having lower annual maximum temperatures. In the year 2021, the total glacial area was below 0.5 km2 in response to warming which showed a maximum annual temperature of more than 55 (F). A positive correlation between the total area and precipitation was observed (Figure 3).
A strong negative correlation (r = −0.59, p < 0.05) between the area of glacier and the annual maximum temperature of Grand Lake 6 SSW station was observed. The areal extent of the glaciers showed a strong positive correlation (r = 0.56, p < 0.05) with precipitation (for each year, from Sep. 1 last year to May. 31 this year) of Grand Lake 6 SSW station. Also, a strong negative correlation (r = −0.46, p < 0.05) to the annual global average temperature from NOAA National Centers for Environmental Information website was observed (Figure 4) (NOAA National Centers, 2021) . Other studies have also found that average temperatures in Rocky Mountain National Park increased by 0.29˚C per decade from 1981 to 2015 (Diaz & Eischeid, 2007; Knowles et al., 2006; Hamlet et al., 2005) .
These correlations show that the retraction of glaciers is a response to the global warming. The precipitation trend in the area is also decreasing, which could also be a reason to the reduction in the total area of the glaciers present in the Colorado Front Range region. Over the past few decades, the climate in the glacial regions of the Northern Hemisphere has become warmer and wetter, especially after 1977. This seems to be an unusual change, possibly man-made (Dyurgerov & Meier, 2000) . In the present study, a significant correlation was found
Figure 3. Line diagram showing variation in the total area of 13 glaciers in the Colorado Front Range region from 1993 to 2021.
between the decreasing areas of cirque glaciers with the decreasing rate of precipitation. Also, a significant correlation was found between the decreasing area of cirque glaciers and the increasing annual maximum temperatures. Thus, glacier retreat might be due to higher summer temperatures since the 1980s (Shahgedanova et al., 2010) . Thus, human-induced climate change could cause cirque glacier retreat, resulting in the disappearance of them around the world.
4. Conclusion
This study fills a gap in the literature by investigating changes in the area for a group of glaciers in the Front Range using the Normalized Difference Snow Index with Landsat imagery from 1993 to 2021. The feasibility of using Landsat
Figure 4. Line diagram showing variation in the total area of 13 glaciers in the Colorado Front Range region from 1993 to 2021 with variation in precipitation, annual regional maximum temperature, and annual global temperature.
images to analyze spatial glacial cover for small cirque glaciers has also been tested, which shows that Landsat images can potentially act as a robust tool to identify the spatial changes of glacial cover. The glacial area calculated using NDSI in the Colorado Front Range region exhibits that the total area of 13 glaciers has decreased significantly from 1993 to 2021. The decreasing areal extent is well correlated with the increasing annual maximum temperature and decreasing rate of precipitation. Thus human-induced increment in global temperature might be an important reason in the decreasing trend of the area of the glaciers in the Colorado Front Range region, USA. We provide a linear regression model, but the non-linear dynamics should be thoroughly explored in the future, which could be used to better explain fluctuations in glacier size due to climatic variables. Besides, potential confounding factors (e.g., land use changes, albedo effect fluctuations) are also needed to be explored in relation to their impact on glacier cover.
Acknowledgements
The authors wish to acknowledge Prof. Dan McGrath, Dr. John Adler, and Prof. Morteza Karimzadeh. Thanks to the geography department in the University of Colorado Boulder, for the help in this study. The authors also thank to the many scientists and practitioners who established, promoted, and improved Google Earth Engine and Google Earth Pro, drove positive scientific and societal outcomes, and the MAXAR for the 41 cm high-resolution image in Google Earth Pro. Finally, we would like to thank the USGS for providing the Landsat imagery.
Data Availability Statement
The data that support the finding of this study are available in Google Earth Engine and Google Earth Pro.