Accelerated retreat of coastal glaciers in the Western Prince William Sound, Alaska

ABSTRACT Analyzing historical maps and Landsat imagery indicates that coastal glaciers in the western Prince William Sound (PWS) have retreated since the end of the Little Ice Age, exhibiting accelerated retreat after the mid-2000s. A multitemporal inventory of 43 glaciers was developed using historical field observations, topographic maps, and Landsat imagery. Area and length measurements are derived from digitized outlines, and center lines are derived from a semi-automatic, geographic information systems (GIS)-based algorithm. Land-based glaciers retreated at a peak rate of 48 m a−1 from the mid-2000s to 2018, more than doubling the average rate of retreat (22 m a−1) for the preceding fifty-year period. From ~1950 to 2018, the total area of land-based glaciers decreased by 228 km2, with 36 percent of the glacier loss occurring after the mid-2000s. Simple upscaling of area and volume changes to unmeasured glaciers across the entire PWS resulted in an estimated aggregate glacier mass loss of 379 Gt, equivalent to a 1.047 mm rise in sea level from the 1950s to 2018. Tidewater glaciers respond asynchronously with differing periods of peak area and length loss and lower average rate of retreat compared to land-based glaciers. Glacier retreat correlates with increased summer and winter temperatures and decreased winter precipitation.


Introduction
Glaciers play an essential role in the study of climate change, and advancing and retreating glaciers provide reliable indicators for climate change at decadal to centennial scales (Adalgeirsdóttir, Echelmeyer, and Harrison 1998;Wiles, Barclay, and Calkin 1999;Calkin, Wiles, and Barclay 2001;Santos et al. 2010;Qiang, Lei, and Chuan-Yan 2014). The recent decline in glacier mass worldwide has important implications for human habitats, including water resource disputes, desertification in areas where glaciers have vanished, and shipping hazards in areas with increased ice loss due to calving from tidewater glaciers (Intergovernmental Panel on Climate Change 2007; Y. Li and Li 2014;Qiang, Lei, and Chuan-Yan 2014). Of particular concern is the relationship between global warming, glacier loss, and sea level rise. Gulf of Alaska glaciers represent an estimated onethird of the global contribution to sea level rise from melting glaciers and ice caps (Arendt et al. 2002;Meier and Dyurgerov 2002). Glaciers in the western Chugach Mountains are among the largest contributors to the rising sea level in the Gulf of Alaska (Berthier et al. 2010).
Comprehensive glacier chronologies and repeat inventories over time for different regions are critical for assessing glacier changes due to the climate or other forces (Nuth et al. 2013;Y. Li and Li 2014;Kienholz et al. 2015). In this study, topographic maps and Landsat imagery were used to document the advance and retreat of land and lake-terminating (land-based) glaciers and coastal tidewater glaciers in the western Prince William Sound (PWS) in the northern Gulf of Alaska ( Figure 1). The forty-three evaluated glaciers terminate on land or in fjords along the PWS and are divided into three subregions: the Sargent Icefield, Spencer-Blackstone Ice Complex, and the western Chugach Mountains. Muir (1902) recorded the earliest observations of glaciers in the area in 1899. Glacier chronologies, extending back to the mid-1700s, exist for some glaciers in the study area (Field 1975;Sturm et al. 1991;Wiles, Barclay, and Calkin 1999;Barclay, Wiles, and Calkin 2003). However, the chronologies only record fluctuations up to the late twentieth century. The aim of this study is to understand how coastal glaciers in the study area fluctuated since the Little Ice Age (LIA) in response to climate change. To fulfill this aim, glacier chronologies, including extent and area, were developed for forty-three coastal glaciers larger than 10 km 2 ( Table 1). The chronologies were based on analysis of topographic maps from ~1950 and Landsat imagery from 1973 and 1975 (1973-1975), 1986, 1994, 2004 and 2006 (2004-2006), and 2018. Further, temporal and spatial changes, including glacier length, area, and mass changes, were analyzed for land-based and tidewater glaciers. Morphometric variables, representing topographic and geometric glacier properties, were also evaluated as potential controls for glacier change. Finally, climate time series and trends were developed and compared with patterns of advance and retreat to understand glacier changes in response to climate change from the 1950s to 2018.

Setting
Glaciers have occupied fjords in the western PWS since the Pleistocene (Wiles, Barclay, and Calkin 1999). The LIA period of glacial advance began across southern Alaska around 1200 AD and continued through the early twentieth century (Porter 1986;Wiles, Barclay, and Calkin 1999). With the exception of some tidewater glaciers, the vast majority of glaciers in Alaska are now retreating (Molnia 2008). Berthier et al. (2018) applied a glacier mass balance approach to detect changes in Alaska's Juneau and Stikine icefields and estimated rapid mass loss from 2000 to 2016. Arendt, Walsh, and Harrison (2009) used repeat airborne altimetry to study forty-six glaciers in Alaska and northwestern Canada and concluded that 75 percent of the glaciers were losing mass at an increasing rate since the mid-1990s, compared to an earlier period beginning in the 1950s to 1970s. The study attributed most glacier mass losses to increased summer (May-September) temperature, with likely contributions from winter (October-April) warming, resulting in increased thermal energy and precipitation as rain instead of snow. Post (1975) described the tidewater glacier cycle as the pattern of slow advance followed by rapid collapse as glaciers pull back from shoals and into deeper water. Meier and Post (1987) studied Columbia Glacier, in the northern PWS, and concluded that "advance and retreat of tidewater glaciers are not related to climate fluctuations in a direct way" (9058) but considered tidewater cycles to be related more to basin and fjord geometry and glacier mass balance. McNabb and Hock (2014) developed a sixty-four-year chronology of glacier length for fifty Alaskan tidewater glaciers to understand the response of tidewater glaciers to climate change. The study relied on digitized glacier outlines derived from topographic maps and Landsat images from 1948 to 2012. The authors concluded that 60 percent of the glaciers are retreating and observed "step change" retreats "where glaciers quickly retreat from a stable position, usually located at a topographic constriction, and stabilize again at another topographic constriction" (McNabb and Hock 2014;165). Tidewater glaciers are influenced by oceanic forces, such as the buoyancy-driven circulation of warm water near the calving front , and respond differently to climate change and often with asynchronous advance or retreat compared to nearby land-terminating glaciers (Post et al. 2011). Brinkerhoff, Truffer, andAschwanden (2017) suggest that ice and sediment dynamics could explain the advance of some tidewater glaciers while most glaciers are retreating.
In the western PWS, the U.S. Geological Survey (USGS) maintains a glacier mass balance monitoring program at Wolverine Glacier, one of the glaciers included in this study. Wolverine Glacier is one of three North American glaciers selected by the USGS for long-term monitoring of glacier response to climate change. Results of long-term mass balance monitoring indicate that Wolverine Glacier is losing mass as a result of warmer summer temperatures (O'Neel et al. 2019). Repeat inventory studies conducted in regions throughout the world show similar trends of recent glacier retreat (Xiao et al. 2007;Błaszczyk, Jania, and Kolondra 2013;Bajracharya, Maharjan, and Shrestha 2014;Osipov and Osipova 2014;Paul and Mölg 2014;Qiang, Lei, and Chuan-Yan 2014;Wei et al. 2014;Winsvold, Andreassen, and Kienholz 2014).
The counterclockwise circulation of the Gulf of Alaska is responsible for the maritime climate in the study area and brings significant moisture into the coastal mountains along the PWS (Beamer et al. 2016). As a result, the region experiences heavy precipitation in the form of rain during the cool summers and snow in the winter. The combination of heavy moisture, cool temperatures, and mountainous topography supports the persistence of coastal glaciers. Mean temperatures for coastal communities in the study area range from −5°C in January to 14°C in July, with average annual snowfall as high as 8,280 mm at Whittier (Alaska Climate Research Center 2019).

Data and methods
Historical topographic maps of the study area provided the spatial data needed to extend glacier length change and area chronologies back to the 1950s. The twenty-one maps obtained for this study are available for download in a georeferenced format from the USGS (2019) allowing for use in geographical information systems (GIS)  (Adalgeirsdóttir, Echelmeyer, and Harrison 1998;Arendt et al. 2002Arendt et al. , 2006McNabb and Hock 2014). Landsat images were accessed from an online service portal (Environmental Systems Research Institute 2019). The images are georeferenced and orthorectified by the USGS, allowing for direct integration into a GIS. The images, at 15-to 60-m resolution, provided the spatial data for the repeat measurement of glacier outlines spanning 1973 to 2018 (Supplementary Table 1). Previous studies provided LIA maximums for eight of the land-based glaciers analyzed in this study (Wiles, Barclay, and Calkin 1999;Barclay, Wiles, and Calkin 2003). The process of developing glacier center lines, necessary for the calculation of glacier length changes, required high-resolution digital elevation data provided by the Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) optical sensor. ASTER Global Digital Elevation Model data (ASTER GDEM.v2) is a product of Japan's Ministry of Economy, Trade, and Industry and the National Aeronautics and Space Administration (30-m resolution, obtained from NASA/METI/AIST/Japan Spacesystems, and U.S./Japan ASTER Science Team 2011).
Annual climate time series were developed for the Cannery Creek, Cordova MK Smith Airport, Seward Airport, and Wolverine Glacier weather stations in the PWS region (Table 2). The climate data are grouped into May-September and October-April averaging periods to represent the summer ablation and winter accumulation seasons in Alaska (Arendt, Walsh, and Harrison 2009;Kienholz et al. 2015). Raw climate data were acquired from the National Oceanic and Atmospheric Administration (NOAA 2019a). Missing data are common in Alaska climate data sets, requiring gap filling in some cases (Supplementary Table 2). Techniques used to fill gaps in monthly records include filling in missing monthly data with daily data averages for the month, when available, and replacing missing months or months with incomplete data (less than twenty-five days of data) with the series average of that month (McAfee, Guentchev, and Eischeid 2013).

Mapping glacier outlines
Randolph Glacier Inventory (RGI Consortium 2017) outline shapefiles, based on imagery from 2006 (Sargent Icefield and Spencer-Blackstone Ice Complex) and 2009 (Chugach Mountain glaciers), provided the initial glacier boundaries. Manual adjustments of the digitized RGI glacier boundaries were based on the interpretation of topographic maps from the 1950s and Landsat images acquired in 1973-1975, 1986, 1994, 2004-2006, and 2018. GIS geoprocessing tools were used to calculate the area for each glacier outline and to measure changes over time (Supplementary Table 3). Glacier outlines are available from the Global Land Ice Measurements from Space (GLIMS) database (Raup et al. 2007; www.glims.org). Issues including the frequency of image collection and persistent cloud cover in the western PWS limited the availability of usable imagery. Several glaciers were excluded from the study due to the presence of debris cover, making it difficult to accuarately digitize glacier boundaries. Each study glacier is identified by a project identification number, including GLIMS and RGI identification numbers as well as glacier name, if available. All topographic maps and Landsat images were either acquired in or projected to UTM Zone 6N with WGS84 datum.

Volume-area scaling
Volume-area power law scaling was used to estimate changes in glacier volume. Bahr, Meier, and Peckham (1997) described the following simple relationship between area and volume: where glacier volume V is estimated from glacier area S, c is a scaling constant in units of length raised to the power (3 − 2γ), and γ is a dimensionless scaling coefficient. Estimates of scaling parameters have been derived empirically from regional data (e.g., Macheret, Cherkasov, and Bobrova 1988;Chen and Ohmura 1990) and theoretical relationships (e.g., Bahr, Meier, and Peckham 1997;Van de Wal and Wild 2001). Among frequently cited volume-area scaling studies, c estimates vary widely by method and region, ranging from c = 0.12 m 3−2γ based on a global glacier model  2001), to c = 0.2055 m 3−2γ from a statistical analysis applied to a large glacier area data set (Radić and Hock 2011). Here, lacking a known c estimate for the study area, the scaling constant c = 0.163 m 3−2γ was derived from the average of the upper and lower ends of range of published c estimates, using γ = 1.375. The scaling equation was applied to the entire population of study glaciers, rather than individual glaciers, to minimize error in volume estimation (Bahr, Pfeffer, and Kaser 2015). Changes in glacier mass, measured in gigatons, were calculated by multiplying ice volumes (km 3 ) by the density of ice and dividing by the time interval. Glacier mass changes and the equivalent rise in sea level (SLE) were calculated for the forty-three measured study glaciers and upscaled to all unmeasured glaciers throughout the PWS, following an approach described by Arendt et al. (2006). The unmeasured glaciers included all 2,828 PWS glaciers available in the GLIMS database (Figure 1), excluding the forty-three measured study glaciers. With the exception of Columbia Glacier, all of the unmeasured glaciers in the PWS were landbased glaciers. Shapefiles of unmeasured glacier outlines, from the mid-to late 2000s, were obtained from the GLIMS database and volume-area scaling was used to calculate an aggregate volume for the entire PWS. In order to estimate aggregate area and volume changes over time, measured glacier area change rates, obtained from analyzing the forty-three study glaciers, were applied to unmeasured glaciers backward to the 1950s and forward to 2018. Columbia Glacier was included in the aggregate volume estimate using a combination of estimates for Columbia Glacier were obtained from the GLIMS database (2007 outline) and manual mapping of Landsat imagery (2018 outline).

Glacier length
Glacier center lines were calculated using a semiautomatic, GIS-based algorithm (Le Bris and Paul 2013), resulting in a center line from the highest point of the glacier to the terminus. Glacier length changes were measured from the intersection of the center line with each glacier terminus, as defined by digitized glacier outlines. Repeat measurements for 1950s topographic maps and the Landsat images acquired in 1973-1975, 1986, 1994, 2004-2006, and 2018, resulted in a glacier length change chronology for each glacier (Supplementary Table 3). For a subset of eight glaciers, length change measurements extended back to the digitized LIA maximum terminus positions identified in previous studies (Table 3).

Morphometric and bathymetric analysis
Following an approach established by Y. Li and Li (2014), six derived morphometric variables representing topographic and geometric glacier properties were used to evaluate potential controls on glacier change. Table 4 includes the list of derived variables. Relative glacier area change (ΔA) was calculated to quantify glacier change from the 1950s to 2018. Glacier area is commonly used to represent change, because GIS allows for easy calculation of glacier area change across time periods (Y. Li and Li 2014). Relationships between ΔA and the six morphometric variables were quantified using multiple regression analysis and interpretation of correlation coefficients (R 2 ) and p values. Meier and Post (1987) observed a relationship between rapid retreat rates and water depth in tidewater glaciers in Alaska. Fjord bathymetric data were obtained from NOAA (NOAA 2019b; https://maps.ngdc.noaa.gov/viewers/bathy metry/) to test this relationship in the study area. Fjord bathymetric data were derived from hydrographic surveys conducted from 1990 to 1997. Based on availability of  depth soundings, mean water depths were calculated for fjord areas between the retreating terminus positions of eight tidewater glaciers, including multiple time periods for Beloit, Coxe, Nellie Juan, Tiger, and Yale glaciers (Table 5).
Regression analysis was applied to assess the relationship between retreat rates and mean water depth and the potential influence of water depth on tidewater glacier fluctuation. Fjord width is another important controlling factor for tidewater glacier retreat and advance (Mercer 1961

Uncertainty
The sources of uncertainty relating to the use of remote sensing in glacier studies are well documented and include satellite image quality and resolution issues, georeferencing issues, manual editing errors (e.g., digitizing outlines), and the analyst's degree of skill in interpreting features such as ice flow divides, snowlines, glacier boundaries, and debris cover (DeBeer and Sharp 2007; Racoviteanu et al. 2009;Paul et al. 2013Paul et al. , 2017McNabb and Hock 2014;Osipov and Osipova 2014;Paul and Mölg 2014;Pfeffer et al. 2014;Qiang, Lei, and Chuan-Yan 2014;Kienholz et al. 2015).
The following study design criteria were applied to minimize potential errors: (1) exclude glaciers with heavy debris cover at margins; (2) avoid satellite images with significant cloud cover or shadows along glacier margins; (3) select late summer (July-September) satellite images to minimize seasonal snow; and (4) obtain georeferenced and orthorectified maps and satellite images, allowing for direct integration into GIS, and check for geolocation errors. Strategies recommended by Paul et al. (2017) were followed to assess the accuracy and precision of glacier outlines and lengths derived from satellite imagery, including the accuracy measures outline overlay, literature values, and multiple digitizing. Outline overlay involves a qualitative analysis of multiple glacier outlines over the source image and helps identify potential sources of error. This method, applied to six study glaciers, revealed differences in outlines in the range of 30-90 m (one to three pixels). The errors were in discrete areas of glacier margins where boundaries between snow and ice, or debris and ice, were difficult to assess.
Literature values provide an additional source of accuracy estimates for digitizing glacier outlines. Uncertainty estimates for manual glacier outline editing range from ~2 percent to 5 percent of the area measured (Paul et  Analyst precision testing in this study followed the methodology recommended by Paul et al. (2017), including the selection of six glaciers of different sizes (10-77 km 2 ) and three iterations of manual outline digitizing for each glacier. The calculated standard deviation ranged from 0.16 to 0.42 km 2 , and the percentage precision ranged from 0.7 to 2.0 percent. For comparison purposes it is important to note that this study avoided glaciers with debris cover along margins, whereas other studies may include all debris-covered glaciers. Bahr (1997) and Bahr, Meier, and Peckham (1997) recognized log-log scaling relationships for glacier volume, area, and length. Establishing area-length scaling relationships allows for the comparison of differences between regions and methods of glacier center line development and provides a simple method for estimating glacier length from area measurements (Le MacHguth and Huss 2014;Kienholz et al. 2015). Area-length scaling for the thirty land-based glaciers in this study resulted in the log-log relationship of L = 1.7836A 0.5416 (R 2 = 0.898), similar to the relationship established by Le , who developed the center line length calculation method applied in this study (Supplementary Figure 1). The large uncertainty resulting from the volume-area power law scaling approach was estimated by calculating glacier volume using the upper (c = 0.2055 m 3−2γ ; Radić and Hock 2011) and lower (c = 0.12 m 3−2γ ; Van de Wal and Wild 2001) ends of the range of published c estimates. The uncertainty associated with the forward and backward extrapolation of area estimates for unmeasured glaciers was calculated as the percentage difference between measured and extrapolated average area change rates for the study glaciers from 1950-1957 to 2004-2006 and 2004-3006 to 2018. Average area change rates were derived from measured study glaciers larger than 10 km 2 . This approach underestimates volume-area scaling results because relative area losses were greater for smaller measured glaciers, and over 98 percent of unmeasured glaciers were smaller than 10 km 2 . The uncertainty was estimated as 3.8 percent for backward (2004-2006 to 1950-1957) and 2.9 percent for forward (2004-2006 to 2018) extrapolation.

Climate analysis
Climate time series were developed for four weather stations in the PWS region to assess the relationship between glacier changes and climate. Time series included temperature and precipitation for each station for annual, ablation season, and accumulation season averages. Detailed climate time series methodology is provided in the Supplementary Material. Previous studies suggest a relationship between abnormally warm summer sea surface temperatures and the onset of rapid glacier retreat (Schild and Hamilton 2013;McNabb and Hock 2014). To further investigate potential factors on tidewater glacier fluctuations, sea surface temperature (SST) data were obtained from NOAA (NOAA 2019c; https://tidesandcurrents.noaa.gov/stations.html) for the Cordova, Seward, and Valdez stations. Mean summer (June-August) SSTs were calculated for Valdez station from 1994 to 2018. Gaps in the water temperature record for Cordova and Seward limited the extent of the mean summer SST analysis to 2015.
Temperature in Alaska is influenced by climate signals including the Pacific Decadal Oscillation (PDO) and shifts in atmospheric circulation in the Arctic. The PDO refers to a cyclical pattern of SST and sea level pressure changes occurring in the North Pacific over periods of ten years or more. North Pacific SST anomalies define the PDO and shift between positive and negative phases. The positive PDO phase is associated with warmer and wetter conditions along the coast of Alaska (Arendt, Walsh, and Harrison 2009). The Arctic Oscillation (AO) refers to the changing strength of the polar vortex (Thompson and Wallace 1998), measured as an index based on sea level pressure anomalies. In the positive AO phase, strong circulation at the North Pole contains the spread of colder air to the midlatitudes, leading to cooler temperatures and increased winter precipitation in Alaska. PDO and AO data from 1950 to 2018 were acquired from NOAA (NOAA 2019d; https://www.ncdc.noaa.gov/teleconnec tions/) and used in regressions to quantify relationships with temperature and precipitation. Time series were developed to provide a visual comparison of glacier changes and significant shifts in PDO and AO index phases.  Table 6). Based on mapped LIA extents (Wiles, Barclay, and Calkin 1999), glacier retreat rates averaged 14 m a −1 from the LIA maximum to the 1950s (Table 3). Average retreat rates for this subset of glaciers peaked after Table 6. Average length and area change metrics for land-based and tidewater glaciers.

Metric
Glacier type Unit 1950s to 1973-1975 1973-1975 to 1986 1986 to 1994 1994 to 2004-2006 2004-2006 (Figures 3b-3c). During this period, 80 percent of land-based glaciers experienced the highest annual area loss rate (Figure 4). The regression between glacier relative area change (1950s-2018) and the six derived morphometric variables, representing topographic and geometric glacier properties, indicates a weak relationship (R 2 = 0.17). Among the six variables, only glacier area showed a significant correlation (p = .03) with glacier relative area change ( Table 7). Mass of the forty-three measured study glaciers declined by 63 ± 10 Gt from the 1950s to 2018, equivalent to a sea level rise of 0.174 mm. Simple upscaling of area and volume changes to unmeasured glaciers across the PWS resulted in an estimated aggregate glacier mass loss for the PWS of 379 ± 38 Gt, equivalent to a 1.047 mm rise in sea level from the 1950s to 2018 (Table 8). The catastrophic retreat of Columbia Glacier accounted for over half (201 ± 11 Gt) of the mass loss for the PWS. Since the 1950s, over 37 percent of the aggregate PWS glacier mass loss occurred after 2004 with an equivalent sea level rise.
Tidewater glaciers behaved asynchronously compared to land-based glaciers, as indicated by differing periods of peak area and length loss and a much lower average rate of retreat. The combined area of study tidewater glaciers declined by 88 km 2 or 6.1 percent since the 1950s (Figure 1; Table 6), with peak losses occurring between 1973-1975 and 1986 (Figure 3b). Tidewater glaciers lost an average of 4.3 percent of total length, compared to a 14.1 percent loss for landbased glaciers since the 1950s (Figure 2). The average retreat rate for tidewater glaciers fluctuated across the study period, with highest rates occurring in 2004-2006 to 2018 and 1973-1975 to 1986 (Figure 3a). The period 1986 to 1994 represents the only study interval with positive average length change rates (2 m a −1 ) for all tidewater glaciers, with four glaciers reaching peak   2003; Figure 6). Nellie Juan Glacier's period of rapid retreat peaked at a rate of 124 m a −1 from 2004-2006 to 2018. During the same period, Barry Glacier retreated over 2 km, into a widening fjord, from a relatively stable position at the head of the inlet shared with Cascade and Coxe glaciers. The glaciers of College Fjord highlight the asynchronous behavior of tidewater glaciers. Since the 1950s, Harvard Glacier advanced slowly, whereas the neighboring Yale Glacier retreated over 5 km (Figure 6).

Climatology
Climate trend analysis included weather station precipitation and temperature data for the four study area weather stations (Supplementary Table 4). Trend analysis results indicate increasing average annual air temperatures in the study area, including an increase from 3.1°C to 4.7°C at Cordova from 1951 to 2017 (Figure 7a). Significant winter warming occurred at PWS maritime stations (+0.9°C) from 1988 to 2017 (Figure 7b). This increase is noteworthy because the average PWS maritime winter temperature is near freezing (−0.5°C). Increasing winter temperatures may expose larger areas of maritime glaciers to temperatures above freezing (Arendt, Walsh, and Harrison 2009). Average summer temperature increased from 5.4°C to 6.1°C from 1968 to 2017 at the higher elevation Wolverine Glacier weather station (Figure 7c). Annual winter precipitation decreased by an average of 319 mm (20 percent decline) at PWS maritime stations from 1976 to 2017 (Figure 7d). Analysis of summer water temperature data obtained from Valdez station indicate a 2.8°C rise from 2008 to 2016 (Figure 8b). Data from Cordova and Seward showed similar trends of increasing summer water temperatures beginning in the late 2000s.
Climate indices provide evidence of significant climate phase shifts in Alaska, with a period of positive PDO from 1977 to 1997 and a large shift to positive AO in 1989 (Arendt, Walsh, and Harrison 2009). Regression analysis was conducted to identify relationships between PDO, AO, air temperature, and precipitation from 1950 to 2018. No significant correlations were identified between PDO and AO. Regression between PDO and annual temperature at Cordova indicates a strong correlation from 1952 to 1989 (R 2 = 0.46; p < .0001); however, the correlation weakens from 1990 to 2017 (R 2 = 0.23; p = .012), suggesting a shift in forces responsible for recent warming in Alaska. Weak correlations were observed between winter PDO and precipitation at Seward (R 2 = 0.12; p = .026) from 1976 to 2017 and Wolverine Glacier (R 2 = 0.15; p = .034) from 1986 to 2015.

Discussion
The majority of coastal glaciers in the study area are retreating, with rapid mass loss occurring from the mid-2000s to 2018. These findings are consistent with reported glacier mass loss in Alaska (Arendt et al. 2006;Molnia 2008;Arendt, Walsh, and Harrison 2009;     Here, measured glaciers in the Sargent Icefield and Spencer-Blackstone Ice subregions lost ~8 percent of their total area from 1986 to 2018. The lower percentage of area loss compared to Yang et al. (2020) is attributed to the exclusion of glaciers smaller than 10 km 2 . Other studies in Alaska (Meier and Dyurgerov 2002), the Northern Alps  (Paul et al. 2004), and the Chinese Tien Shan (K. Li et al. 2011;Y. Li and Li 2014) suggest that glacier area as a key indicator of glacier loss, with smaller glaciers tending to have greater relative reductions. Arendt et al. (2006) used airborne altimetry measurements and multiple extrapolation methods, including special considerations for tidewater glaciers, to estimate a regional contribution of 0.02 mm yr −1 SLE and a net balance rate of −7.4 ± 1.1 km 3 yr −1 water equivalent for glaciers covering 9,300 km 2 of the western Chugach Mountains and Talkeetna Mountains from 1950-1957 to 2001-2004. Here, simple upscaling to a similar total aggregate glacier area of 9,600 km 2 , including glaciers in the western Chugach Mountains, Sargent Icefield, and Spencer-Blackstone Ice Complex, results in an estimated rate of 0.013 mm yr −1 SLE and a net balance rate of −4.6 ± 0.5 km 3 yr −1 water equivalent from 1950-1957 to 2004-2010. The asynchronous glacier cycles observed among tidewater suggest the influence of factors other than climate, such as fjord geometry, water depth, mass balance, and sediment dynamics (Meier and Post 1987;Brinkerhoff, Truffer, and Aschwanden 2017). Warming waters in the PWS may have contributed to the increased average rate of tidewater glacier retreat observed since the mid-2000s (Figures 8b, 8e), and the rapid collapse of Nellie Juan and Barry glaciers. Regionally, SST anomalies throughout the Gulf of Alaska increased by up to 3°C to 4°C since 2013 (Campbell 2018). Other studies from Alaska and Greenland stress the importance of ocean water temperatures as a controlling factor for submarine glacial melting and calving rates (Motyka et al. 2003Post et al. 2011).
Earlier studies linked rapid retreat of glaciers in the Gulf of Alaska to water depth and emphasized the importance of nonclimatic controls on tidewater glacier fluctuations (Post 1975;Meier and Post 1987). Here, regression analysis between tidewater glacier retreat rate and mean water depth (Figure 9a) results in a significant positive relationship (R 2 = 0.58; p < .001). Water depth may explain the rapid retreat of Yale Glacier, occurring from 1973 to 1986. During this period, Yale Glacier retreated over the deepest waters in Yale Arm at a rate of 179 m a −1 (Table 5). Further investigation into fjord geometry suggests that fjord width may be another controlling factor for tidewater glacier advance and retreat. From the mid-2000s to 2018, all four advancing tidewater glaciers moved into a narrowing fjord, whereas Nellie Juan and Barry glaciers experienced the most rapid retreat in fjords that were rapidly widening (Figure 9b). The relationships presented here suggest that factors other than climate are influencing tidewater glacier fluctuations. McNabb and Hock (2014) determined that tidewater glacier advance and retreat patterns in Alaska are highly variable and noted a lack of evidence for a coherent regional-scale retreat. Here, the lower average rates of tidewater glacier retreat and area loss compared to landbased glaciers are attributed to nonclimatic controls and the variability in tidewater glacier cycle position for the thirteen studied tidewater glaciers.
Increasing temperatures and decreasing winter precipitation explain overall glacier recession in the study area. The accelerated retreat after 2004-2006 is likely attributable to an unfavorable combination of temperature and precipitation changes that occurred before and after [2004][2005][2006]. Figure 8 shows the relationship between the average length change rates, temperature precipitation trends for summer and winter seasons, summer SSTs, and climate indices. Summer temperatures and SSTs increased from ~2012 to 2017. To further drive accelerated melting after 2004-2006, winter temperature increased sharply at all stations from ~2012 to 2016, and summer precipitation increased at all stations from ~2009 and ~2013. Preceding this period of increased summer precipitation, warming, and accelerated glacier retreat, a steep decline in winter season precipitation began around 2002 and continued to ~2010 (Cannery Creek and Seward).
The brief shifts to negative PDO and positive AO after 1988 coincide with a period of positive average length change rates for all tidewater glaciers (Figures  8c-8e). However, this relationship does not extend to the recent period of increased temperature and accelerated glacier retreat. After 2004, PDO was strongly negative with the exception of a short period of positive PDO from 2014 to 2016 (Figure 8c). Arendt, Walsh, and Harrison (2009) suggested that forces other than climate indices may be responsible for warming after 1999, when AO began to trend neutral and PDO shifted to the negative phase.

Conclusions
Analysis of topographic maps from the 1950s and Landsat imagery from 1973-1975, 1986, 1994, 2004-2006, and 2018 supported the development of repeat glacier inventories for forty-three glaciers in the western Prince William Sound, Alaska. Coastal glaciers in study area have been retreating since the end of the LIA, with a period of accelerated retreat and mass loss from the mid-2000s to 2018. Simple upscaling to the entire PWS results in an estimated aggregate glacier mass loss of 379 Gt, equivalent to a 1.047 mm rise in sea level from the 1950s to 2018. Glacier retreat correlates with climate trends including increased summer and winter temperatures and decreased winter precipitation over the last several decades. Sequencing of declining winter precipitation before 2004-2006, followed by increasing summer precipitation and summer and winter temperatures, may have forced accelerated melting after the mid-2000s. The lack of significant statistical relationships between topographic and geometric glacier properties and relative glacier area change further supports the suggestion that climate is the key driver for land-based glacier loss in the study area.
Tidewater glaciers behaved asynchronously compared to land-based glaciers, as indicated by differing periods of peak area and length loss and a much lower average rate of retreat (Figure 3). However, though some tidewater glaciers are advancing, the combined area of all tidewater glaciers in the study declined by approximately 6 percent since the 1950s. This overall decline in tidewater glacier area is likely driven by factors other than climate, such as fjord geometry, bathymetry, and warming summer SSTs. However, climate may have triggered cycles of retreat and indirectly contributed to the overall loss in tidewater glacier mass because warming waters in the PWS likely increased calving rates and the rapid collapse of Barry, Nellie Juan, and Yale glaciers. The tidewater glacier cycle and high variability in advance and retreat patterns for the thirteen study tidewater glaciers explains the lower average rates of tidewater glacier retreat and area loss compared to land-based glaciers.
Expanding the study to the eastern PWS and along the Gulf of Alaska may provide insights regarding the extent of recent accelerated retreat documented in this study. Incorporating remote sensing-based or physical measurements of glacier thinning and volume changes would provide a more accurate estimate of glacier response to climate change in the study area. Further glacier-specific studies are warranted to better understand the complex interaction of climate, fjord geometry, and tidewater glacier mass balance.