Seven Decades of Coastal Change at Barter Island, Alaska: Exploring the Importance of Waves and Temperature on Erosion of Coastal Permafrost Bluffs

: Observational and the number of niche-forming and block-collapsing episodes associated with higher air and water temperature, more frequent storms, and longer ice-free conditions in the Beaufort Sea.


Introduction
Processes driving coastal change in permafrost environments can be fundamentally reduced to mechanical and thermal processes [1,2], and shorelines anchored by permafrostrich substrates and those influenced by the presence of sea ice, are highly vulnerable to the effects of climate change. Observed and projected increases in the duration of the ice-free season, warming of land, air, and ocean temperatures, changes in wave intensity Barter Island is a remnant tundra island that sits between two large deltas of the Jago and Okpilak-Hulahula Rivers. It is one of the largest islands along Alaska's coast, with a total area of nearly 1600 hectares (six square miles). The geological history of the island is Barter Island is a remnant tundra island that sits between two large deltas of the Jago and Okpilak-Hulahula Rivers. It is one of the largest islands along Alaska's coast, with a to-Remote Sens. 2021, 13, 4420 4 of 25 tal area of nearly 1600 hectares (six square miles). The geological history of the island is not well known, but surficial coastal plain sediments in the area are thought to be composed of reworked eolian, alluvial, fluvial, and marine deposits of the Quaternary Gubik Formation, covered by a 2-3 m thick layer of late Pleistocene and Holocene peaty/muddy thaw-lake sediment [11][12][13]. Barter Island is anomalously high compared to the surrounding coastal features, reaching up to 16 m high in the central part of the island ( Figure 1B). Broad low-lying (<2 m high) sand and gravel spits extend east and west from the topographically higher tundra hinterland ( Figure 1C). The community of Kaktovik and adjacent U.S. Air Force radar site sit on the north shore of the island near the eastern end of a nearly 3 km stretch of eroding coastal permafrost bluffs. The bluffs range in height from a few meters to more than ten meters high and consist of a complex sequence of poorly consolidated material ranging from dense clay, interbedded and massive sand and gravel, massive ice, which has been interpreted as buried glacial ice [14][15][16], wedge ice, segregated ice layers, massive ice, and thermokarst cave ice, all capped with a surface peat layer [13,16,17]. Barter Island, like the entire North Slope of Alaska, is in a zone of continuous permafrost, with a shallow active layer (historically generally less than 1 m deep) that is subject to thawing during the summer season. Narrow (<20 m wide) beaches front most of the coastal bluffs. Bernard Spit, a low-lying barrier island to the northeast, protects much of Barter Island's eastern spit from incident wave energy. The island is subject to overwash and breaching during periods of high waves and elevated water level. There is an apparent absence of ice-bonded permafrost in Kaktovik Lagoon [18].
The Beaufort Sea is typically covered with sea ice from approximately October through June. Prevailing winds are easterly and average near 5.8 m/s with little annual variation [19,20]). The region is microtidal with a mean diurnal range of 18 cm [21]; however, water levels can become elevated or depressed up to several meters due to winds and low-pressure systems; westerly winds tend to elevate water levels, whereas easterly winds tend to lower water levels [4, 22,23]. Wave directions are predominantly from the east, but the largest waves are associated with the passage of large low-pressure systems that generally approach from the west. Studies of projected climate patterns indicate a poleward shift in storm tracks [24] and particularly along the Barter Island coast [4]. The near E-W trending Barter Island shoreline coincides with a nearly balanced E-W longshore transport regime that drives sediment redistribution along the open coast. This is reflected by the large sand spits that extend both east and west of the island.
Coastal geomorphic processes are, for the most part, limited to three ice-free months of the year between late July and the onset of shorefast ice in October. Starting in late July to mid-August, low-pressure cells originating in the Pacific Ocean or Chukchi Sea traverse eastward across Alaska and the Beaufort Sea bringing with them high winds that can generate storm surges possibly as high as 3 m and waves as high as 6 m [4, 12,[25][26][27]. The magnitude of littoral transport on the beaches and barrier islands in this area is thought to be low along this coast compared to lower latitude coasts due to relatively low wave energy and short open water season [12]. Sediment transport directions are bimodal with dominant easterly waves generating limited longshore sediment transport from east to west and infrequent but larger westerly waves generating longshore sediment transport from west to east.
Similar to permafrost bluffs elsewhere in the Arctic, primary bluff failure modes are a combination of thermal denudation on the bluff face throughout the summer and thermomechanical removal of debris material fronting the bluffs and niching of the base associated with waves and high-water levels (wind or storm surge) followed by rotational slumping and block collapse [1,2,[28][29][30][31]. Ground ice typically occupies 60-90% of the volume of near-surface deposits [16,32], and it is a major factor contributing to the high coastal erosion rates [33]. Retreat rates appear to be largely dependent on ice content, the frequency and intensity of storms, runup elevation, and seawater and air temperatures [34][35][36]. The pattern of change is predominantly landward retreat of the top of the bluffs, removal of the debris apron and subsequent niching at the base of the bluffs, followed by continued erosion of the bluff face and deposition of debris at the base of the bluff [29].

Methods
Methodologies used to quantify rates and patterns of landscape change and three potential physical drivers of change are outlined in Figure 2. Shoreline and bluff edge positions were delineated from a combination of historical maps, aerial photography, declassified satellite photography, and very-high resolution satellite imagery spanning a 73 year period between 1947 and 2020. Spatial resolution of the datasets ranges from 0.08 to 4.1 m (Table 1). Shoreline data comprise 21 individual years spanning 1947 to 2020; however, there are only 6 years with complete shoreline coverage over the study area. Bluff edge data comprise 22 individual years ranging from 1950 to 2020. A useable bluff edge position could not be determined from the 1947 nor 1987 T-sheet data. Images were georeferenced and rectified using map coordinates, ground survey points, and/or image-to-image correlation of photo-identifiable ground features. A total root mean square rectification error was determined for each image based on the residuals. The top edge of the bluff and the instantaneous land-water shoreline were manually digitized in a GIS, and a total position uncertainty of the feature was determined as the quadrature sum of all uncertainties, including image resolution, georeferencing/rectification, and feature digitization [9] (Table 1).

Shoreline and bluff
Bold text indicates shoreline datasets used in this analysis; grayed text indicates bluff and shoreline positions are included in [29] but not used in this analysis.
Change rates were calculated every 10 m alongshore using the Digital Shoreline Analysis System 5.0 [37] using either the linear regression or endpoint rate calculation methods, depending on the number of features selected for analysis (greater than 2 or equal to 2, respectively). Bluff change rates were calculated for 4 multidecadal time periods: 1950 to 2020, 1950 to 1975, 1975 to 2000, 2000 to 2020, and 21 semiannual to annual time periods (Table 2). Regionally averaged rate uncertainties were calculated using the reduced n approach on long-term linear regression rates [38,39] and the mean of the endpoint rate uncertainties (EPRunc) on individual transects as described in [37]. Shoreline change rates were calculated along the open coast from the east end of Bernard Spit to the west extent of Barter Island for the 6 years with data extending across the entire study area (1947,1979,2004,2009,2018, and 2020; Figure 1, Table 1). Additional details on methodology and uncertainty estimates can be found in [9,37]. Bluff edge and shoreline vectors, rate calculation data, and error and uncertainty statistics at individual transects are available at [40].
Attribution of bluff recession to three physical drivers was evaluated by computing cumulative wave power, ocean surface temperatures, and ambient air temperatures, and comparing these to bluff recession rates. Whereas process-based models are preferred for their granularity and better accuracy, here we employed wave power (propagated to the inshore region) and sea surface temperatures as indicators of thermoabrasion, and ambient air temperatures as proxies for thermodenudation.
Wave power (WP) was computed with, where ρ is the density of water (=1028 kg/m 3 ), g is the gravitational acceleration (=9.83 m/s 2 at N70.1 • ), H s is the significant wave height (units of meters), and T e is the energy period (units of seconds), herein estimated with 0.538·T 01 , based on the JONSWAP wave spectrum [41]. Akin to the concept of 'growing-degree-days' frequently used in agriculture, whereby differing threshold temperatures are used to determine accumulated heat units for different crops, we compute 'positive-degree-days' (PDD) as a possible indicator of thermodenudation. The PDD was computed as the cumulative sum of daily mean air temperatures above 0 • C. Sea surface temperatures were similarly treated to compute PDSST. For evaluation of changes since 1979, WP, PDD, and PDSST were computed annually, whereas for evaluation with respect to bluff recession rates, WP was cumulated between subsequent timepoints of bluff position measurements, whilst PDD and PDSST were accumulated between the most recent of either the preceding measurement date or start of the year (i.e., thermal denudation was 'reset' to the current measurement date if multiple years passed between measurements).
All data used to estimate WP, PDD, PDSST, and duration of open water season were extracted from the ERA5 reanalysis (ERA5 | ECMWF). ERA5 is the 5th generation European Centre for Medium-Range Weather Forecasts (ECMWF) reanalysis of the global climate. It was developed by combining model data with vast observations from across the world into a globally complete and consistent dataset. The reanalysis dataset was created using a 4D-VAR sophisticated data assimilation method as part of the new ECMWF Integrated Forecasting System (IFS) CY41R2; it simulates hourly atmospheric fields at 0.25 • (~30 km) spatial resolution from 1979 onwards. The reanalysis additionally includes ocean wave parameters (at a 0.50 • (~56 km) resolution) derived from a fully coupled atmosphere-oceanwave model (ecWAM) that incorporates assimilation of satellite radar altimeter-derived wave height data from 1991 to present. Air and sea surface temperatures employed in this study were extracted from the grid cell closest to Barter Island. Wave parameters (heights, periods, and incident directions) were extracted from a grid cell approximately 30 km offshore in~20 m water depth. Because wave energy is a dominant driver of coastal erosion, and because ERA5, though relatively fine scale considering the global extent of the reanalysis product, is not able to resolve variations of wave energy within the inner nearshore region, the ERA5 wave time series was downscaled to the inner nearshore region using a local-scale wave model. pean Centre for Medium-Range Weather Forecasts (ECMWF) reanalysis of the global climate. It was developed by combining model data with vast observations from across the world into a globally complete and consistent dataset. The reanalysis dataset was created using a 4D-VAR sophisticated data assimilation method as part of the new ECMWF Integrated Forecasting System (IFS) CY41R2; it simulates hourly atmospheric fields at 0.25° (~30 km) spatial resolution from 1979 onwards. The reanalysis additionally includes ocean wave parameters (at a 0.50° (~56 km) resolution) derived from a fully coupled atmosphere-ocean-wave model (ecWAM) that incorporates assimilation of satellite radar altimeter-derived wave height data from 1991 to present. Air and sea surface temperatures employed in this study were extracted from the grid cell closest to Barter Island. Wave parameters (heights, periods, and incident directions) were extracted from a grid cell approximately 30 km offshore in ~20 m water depth. Because wave energy is a dominant driver of coastal erosion, and because ERA5, though relatively fine scale considering the global extent of the reanalysis product, is not able to resolve variations of wave energy within the inner nearshore region, the ERA5 wave time series was downscaled to the inner nearshore region using a local-scale wave model. To propagate WP to the inner nearshore (here defined to be within the 5 m isobath and shoreline) a SWAN wave model consisting of a curvilinear grid with nominal resolution of ~200 m was built. Bathymetry was populated with data from single-beam surveys collected in 2011 [4] and IBCAO [42]. Although highly dynamic, the position of Bernard Spit was kept constant with a footprint that was representative of the 2010 decade. The model was run for 2000 separate sea states representing possible combinations of winds and waves at the ERA5 offshore point for the 1979-2019 time period. The representative To propagate WP to the inner nearshore (here defined to be within the 5 m isobath and shoreline) a SWAN wave model consisting of a curvilinear grid with nominal resolution of 200 m was built. Bathymetry was populated with data from single-beam surveys collected in 2011 [4] and IBCAO [42]. Although highly dynamic, the position of Bernard Spit was kept constant with a footprint that was representative of the 2010 decade. The model was run for 2000 separate sea states representing possible combinations of winds and waves at the ERA5 offshore point for the 1979-2019 time period. The representative sea states were determined using a clustering technique developed by [43]. WP at five locations fronting Remote Sens. 2021, 13, 4420 9 of 25 the bluffed Barter Island coast was computed from all 2000 SWAN runs and subsequently divided by WP at the offshore input point to derive WP reduction factors, α, where the i subscript denotes individual sea states, OF is the offshore (ERA5) point, r is one of the five inshore points along the 2 m isobath fronting one of the five bluff sections ( Figure 3A), N is the total number of sea states, and mean values are denoted by the overbar. The mean valued reduction factors were applied to the hourly ERA5 time series of computed WP to generate inshore WP between subsequent measurements of bluff top positions (denoted with a superscript T).

Shoreline Change
The sand and gravel shorelines of Barter Island and adjacent barrier island, Bernard Spit, showed a consistent pattern of southerly (landward) retreat and westerly extension through time ( Figure 3A). Between 1947 and 2020 the east and west ends of Bernard Spit eroded and/or migrated to the south while the central portion of the island remained relatively stable. The western terminus of the island retreated up to 470 m and extended nearly 980 m to the west. The total area lost (P1 in Figure 3B) and gained (P2 in Figure 3B) on the western terminus were nearly equivalent (332,000 and 329,000 m 2 , respectively), indicating a relative balance of subaerial sediment redistribution. A volumetric comparison of the observed change was not possible because of lack of available elevation data prior to 2009. The eastern terminus of the island was breached, and a persistent inlet formed, sometime between 1979 and 2004, ultimately resulting in over 600 m of southerly migration of the eastern terminus by 2020.
On Barter Island, maximum beach widths were present at the flanks of the island where east-west trending spits formed, by way of longshore transport driven by variations in wave direction at the terminus of the spits. The presence and westward migration of Bernard Spit has increasingly provided shelter from incident wave energy to eastern Barter Island, resulting in relative stability of the eastern spit over the study period. The beach fronting the community of Kaktovik and adjacent lagoon, began prograding northward, and progressively westward, beginning sometime between 1987 and 2000. By 2020 a relatively persistent beach, over 100 m wide, was present as far west as the radar station ( Figure 3A). The shoreline fronting the coastal bluffs to the west of the radar station generally retreated at a rate similar to the eroding bluffs, while maintaining a narrow, but seasonally variable width of between 0 and approximately 30 m [29]. Similar to Bernard Spit, Barter Island's western spit, changed considerably over the study periods, retreating southward as much as 277 m, eroding into and depositing on top of the low-lying vegetated wet-sedge habitat, and extending over 1 km to the southwest ( Figure 3B).
The rate of migration of the western termini of the western Bernard and Barter Island spits were estimated by measuring the linear distance between the westernmost point on the spit in 1947 and the westernmost point in each subsequent year ( Figure 4). Results of this approach showed relatively similar patterns and rates of change on both spits but with a fundamental shift in behavior between 1987 and 2000. Prior to 1987 rates of westward migration of the west Barter Island spit were slightly higher than west Bernard Spit, with rates of 7 and 9 m/yr, respectively. Between 1987 and 2000, the western node of the Barter Spit moved to the north east by approximately 400 m before continuing a near linear and continuous westward migration of near 38 m/yr. The rate of migration rate similarly increased substantially on west Bernard Spit to 32 m/yr but without a notable shift in the node of migration.

Coastal Permafrost Bluff Change
Barter Island's ocean exposed coastal permafrost bluffs retreated continuously over the study period, reaching a maximum distance of 163 m, and average distance of 114 m, between 1950 and 2020 ( Figure 5). The mean long-term bluff retreat rate of 1.6 ± 0.1 m/yr (range 0.3 and 2.2 m/yr) was similar to the long-term mean shoreline retreat rate measured along the Alaskan Beaufort Sea Coast (1.8 ± 0.02 m/yr, 1947-2012) [9] and nearly two times higher compared to mean shoreline retreat rate on the adjacent Canadian Yukon Beaufort Sea coast (0.7 ± 0.2 m/yr; 1951-2011) [44]. The long-term record at Barter Island was also punctuated by individual years with mean annual retreat rates that were more than 4 times higher than the long-term rate ( Table 2). The highest rate of retreat at a single transect (21.3 ± 1.9 m/yr) was measured between July 2012 and August 2013 with a mean rate of retreat for all transects measured during the same time period of 6.6 ± 1.9 m/yr. The region of maximum retreat over the long term was measured near the center of the study area, from approximately 1060 to 1840 m from the eastern boundary of the bluffs. Minimum long-term retreat distances were measured along a ~350 m long section of bluffs on the eastern end of the study area fronted by a prograding beach where retreat has been

Coastal Permafrost Bluff Change
Barter Island's ocean exposed coastal permafrost bluffs retreated continuously over the study period, reaching a maximum distance of 163 m, and average distance of 114 m, between 1950 and 2020 ( Figure 5). The mean long-term bluff retreat rate of 1.6 ± 0.1 m/yr (range 0.3 and 2.2 m/yr) was similar to the long-term mean shoreline retreat rate measured along the Alaskan Beaufort Sea Coast (1.8 ± 0.02 m/yr, 1947-2012) [9] and nearly two times higher compared to mean shoreline retreat rate on the adjacent Canadian Yukon Beaufort Sea coast (0.7 ± 0.2 m/yr; 1951-2011) [44]. The long-term record at Barter Island was also punctuated by individual years with mean annual retreat rates that were more than 4 times higher than the long-term rate ( Table 2). The highest rate of retreat at a single transect (21.3 ± 1.9 m/yr) was measured between July 2012 and August 2013 with a mean rate of retreat for all transects measured during the same time period of 6.6 ± 1.9 m/yr. The region of maximum retreat over the long term was measured near the center of the study area, from approximately 1060 to 1840 m from the eastern boundary of the bluffs. Minimum long-term retreat distances were measured along a~350 m long section of bluffs on the eastern end of the study area fronted by a prograding beach where retreat has been

Long-Term Decadal Scale Bluff Change
Comparison of bluff change rates calculated over near similar length decadal time periods, constrained by limited data prior to 2000 (1950-1975, 1975-2000, and 2000-2020), showed a steady increase in retreat rates through time, with a near 3-fold increase in mean (and maximum) rates of retreat over the last two decades compared to the two previous 25 year time periods ( Table 2)

Long-Term Decadal Scale Bluff Change
Comparison of bluff change rates calculated over near similar length decadal time periods, constrained by limited data prior to 2000 (1950-1975, 1975-2000, and 2000-2020), showed a steady increase in retreat rates through time, with a near 3-fold increase in mean (and maximum) rates of retreat over the last two decades compared to the two previous 25 year time periods ( Table 2). The mean 1950 to 1975 retreat rate of 1.1 ± 0.2 m/yr increased slightly, but significantly, to 1.4 ± 0.2 m/yr from 1975 to 2000 (similar to the 70 year mean of 1.6 ± 0.1 m/yr), then nearly doubled to a rate of 2.7 ± 0.1 m/yr from 2000 to 2020 ( Figure 5). Variability also increased considerably during this final time period, as indicated by an increase in standard deviation of the rates from 0.4 to 1.4. Figure 5C shows the spatial variability in retreat rates across the study area over these same decadal time periods and compared to the 70 year long-term rate. Change rates from the earliest time period, 1950-1975, were consistently lower than the long-term rate. The construction of the landfill was responsible for the accretion around 750 m, and the highest rates of retreat were measured around 1500 m from the east end of the bluffs. Change rates from 1975 to 2000 were similar to the 70 year mean rate, and the highest rates of retreat were measured around 1300 m from the east end of the bluffs. Between 2000 and 2020, change rates were comparatively much higher across the study area except for the easternmost bluffs where retreat rates approached zero. The highest rates of retreat were measured approximately 2000 m from the east end of the bluffs. Similar to the western termini of Bernard Spit and Barter Island spit, the apex or node of the island, which was generally the region of highest retreat rates, also migrated to the west, approximately 800 m since 1950, resulting in an overall straightening of the shoreline line, possibly caused by divergence in longshore transport due to changes in regional or local incident wave directions (Figure 4).
Quantifying the rate or pattern of this apparent migration in the node of maximum erosion was not as straightforward for the coastal bluffs as compared to the sand and gravel shorelines, likely because of the more complex relationship between coastal geology and physical processes driving erosion. The distance between the location of maximum retreat measured over each interannual time period and the location of overall maximum bluff retreat measured over the study period (1950 to 2020) is plotted in Figure 6. No significant temporal trends were apparent over the time series; however, a weak but relatively consistent trend, similar to the spits, was apparent after 2010 through 2020 (Figure 4).
A quartile distribution of total (70 year) bluff retreat distance at individual transects illustrates the spatial variability across the study region as distinct regions with similar magnitudes of change: Central bluffs (Q4; 136 to 163 m), Western bluffs (Q3; 107 to 136 m), and Eastern bluffs, which are further subdivided by moderate (Q3; 102 to 136 m), low (Q2; 59 to 102 m), and very low (Q1; 14 to 59 m) retreat rates ( Figure 7A). The cumulative distance of change for all transects at each bluff-edge observation date illustrates the temporal variation of the bluff over the 70 year study period ( Figure 7B). Fitting linear regression lines to the full dataset, the quantile divisions, and Transect 136, the transect with the highest measured retreat rate from 2019 to 2020, illustrates how each of these spatially varying populations are retreating compared to the average rate of the entire bluff. Applying a 2nd order polynomial linear regression improved the fit of the curves by 0.2 to 4%, particularly for Q3 and Q4 locations, supporting other observations that retreat rates are accelerating through time.

Short-Term Interannual Bluff Change
Short-term interannual bluff change rates showed large temporal and spatial variations, which are characteristic of highly dynamic coastal environments. Although comparison of endpoint rates of change calculated on temporally adjacent bluff positions can provide useful information on interannual variability of the natural system, care must be taken when extrapolating these observations in the context of interpreting meaningful or sustained change to coastal behavior or response to environmental drivers. In the dynamic coastal environment, rates of change calculated over a shorter period of time can be highly skewed by annual rates that are anomalously high (or low) but will generally always be higher in magnitude and variability (higher standard deviation) compared to longer-term measurements, due to the averaging effect over time. Here, we present rates and patterns of interannual change, and in the context of longer-term averages, to highlight significant changes in temporal and spatial trends.     Annual to semiannual (hereafter, interannual) bluff change rates were calculated between adjacent years of bluff position data for a total of 21 comparison pairs, ranging between 0.8 and 21.3 years, with a higher frequency after the year 2000 and near annual comparisons after 2006 (Table 2, Figure 8). Prior to the year 2000, mean short-term bluff retreat rates for four of the five time periods were lower (0.7 ± 1.1 to 1.5 ± 0.2 m/yr) but not significantly different than the long-term 70 year mean retreat rate of 1.6 ± 0.1 m/yr. The exception was the earliest time period, 1950-1955, which was higher (1.9 ± 0.9 m/yr) but also not significantly different than the mean long-term rate. After the year 2000, mean short-term rates ranged from 0.9 ± 2.3 to 6.6 ± 1.9 m/yr, and 11 of 16 years had rates higher than the 70 year mean. After 2015, mean annual change rates were from 1.4 to 2.5 times greater than the 70 year mean (Table 2, Figure 8).
Change rates were calculated for five, 4 year time-averaged periods between 2000 and 2020 to reduce the high variability in annual change rates inherent in the post-2006 data. Mean bluff retreat rates for these 4 year time periods ranged from 1.8 to 3.5 m/yr, were from 1.1 to 2.2 times higher than the 70 year mean retreat rate and showed a consistent and increasing trend through time where retreat rates doubled over the 20 year period ( Figure 8). These results suggest that a steady increase in bluff retreat rates observed since 2000, and consistently since 2015, indicate an acceleration in bluff erosion that is independent of large spatial and temporal variations observed on an annual basis.
Rates of retreat measured on individual transects across the study area for each interannual time period ( Figure 6) illustrate the extremely high spatiotemporal variability of the bluff change dataset. Interannual maximum retreat rates, shown by the black circles in Figure 6  . Mean annual retreat rate and standard deviation for the time periods between adjacent bluff datasets (gray), with mean long-term decadal rates (black) and 4 year mean rates from 2000 to 2020 (red) show highly variable, but consistent increase in retreat rates over the study period.

Assessment of Physical Drivers: Waves and Temperatures
Cumulative wave power (WP) calculated at the western bluffs site, the duration of the open-water season (OWS; within at least 70 km of the coast) derived from the ERA5 reanalysis, and positive degree days and SSTs (PDD and PDSST) derived from the ERA-I reanalysis that were normalized to the 30 year mean (1979-2009) highlight long-term variability in these potential drivers of coastal change ( Figure 9A,B). WP and OWS tracked similarly over the record, illustrating their interdependence with one another [3,45]. The two showed a consistent negative anomaly, or phase, prior to 1987, followed by a mixed phase with both positive and negative excursions from the 30 year mean. A positive phase began in 2003 that was continuous and generally increasing through time, with the excep- Figure 8. Mean annual retreat rate and standard deviation for the time periods between adjacent bluff datasets (gray), with mean long-term decadal rates (black) and 4 year mean rates from 2000 to 2020 (red) show highly variable, but consistent increase in retreat rates over the study period.

Assessment of Physical Drivers: Waves and Temperatures
Cumulative wave power (WP) calculated at the western bluffs site, the duration of the open-water season (OWS; within at least 70 km of the coast) derived from the ERA5 reanalysis, and positive degree days and SSTs (PDD and PDSST) derived from the ERA-I reanalysis that were normalized to the 30 year mean (1979-2009) highlight long-term variability in these potential drivers of coastal change ( Figure 9A,B). WP and OWS tracked similarly over the record, illustrating their interdependence with one another [3,45]. The two showed a consistent negative anomaly, or phase, prior to 1987, followed by a mixed phase with both positive and negative excursions from the 30 year mean. A positive phase began in 2003 that was continuous and generally increasing through time, with the exception of years 2013, 2019, and 2020 ( Figure 9A). PDD and PDSST generally tracked together over the analysis period and, similar to WP and OWS, showed a negative phase prior to 1987, a mixed phase that extended through 2009, followed by a sharp increase in temperature that persisted from 2010 to 2018 ( Figure 9B).
Despite a large year-to-year variability in the first and, to a lesser extent, the last day of open water, the overall duration of the OWS near Barter Island has more than tripled from~40 days in the early 1980s to~140 days by 2020 ( Figure 9C Computed ratios of inshore to offshore wave power along the five distinct bluff sections ( Figure 7A) are shown in Figure 10. The wave energy reduction coefficients were computed separately for five offshore incident wave direction bins (ranging from east (90) to west (270)) and show that west and northwesterly incident waves, once propagated to the shore (~2 m water depth), underwent greater loss of energy along the easternmost bluff section compared to the other bluff sections. For these wave conditions, inshore wave energy at the easternmost bluff section was~60% of the offshore energy, whereas along the western section, the wave energy was~80% of the offshore energy (a reduction of 20% compared to 40% for the easternmost section). For the most easterly incident waves, the presence of Bernard Spit blocks much of the wave energy from reaching any section of the bluffed Barter Island coast. The blocking effect decreases with increasingly westerly incident waves as shown with the relative increase in α(alpha) along the central and western bluffs. Whereas Bernard Spit does not directly block north to west incident wave energy, the presumed longshore sediment transport and westward migration of Barnard Spit's west terminus has led to a wider and milder sloping seabed fronting the eastern bluffs, resulting in greater wave refraction and energy dissipation along that section of coast compared to the central and western sections (see bathymetric contours in Figure 3B).
Median wave power values were applied to a continuous hourly hindcast (1979 through 2019) of offshore ERA5 wave conditions to estimate annual cumulative wave power fronting each of the five bluffed coastal sections ( Figure 11). The presence of landfast and pack ice during most of the year limited the analysis period to late May through early December. Annual variations in timing of spring season sea ice breakup in conjunction with wind energy of sufficient magnitude to generate waves can be inferred from the shift in the gray lines along the x-axes of Figure 11A The results indicate that inshore wave power has increased four-fold over the past decade compared to the long-term mean, and that the increase is continuing (purple line > blue line), commensurate with the increasing bluff retreat rates. WP and OWS data calculated from daily means derived from the ERA5 reanalysis data, with the wave power downscaled to the 'west bluff' nearshore point location using the inshore:offshore ratios presented in Figure 10. PDD and PDSST data were calculated from daily means derived from the ERA-I reanalysis data.  WP and OWS data calculated from daily means derived from the ERA5 reanalysis data, with the wave power downscaled to the 'west bluff' nearshore point location using the inshore:offshore ratios presented in Figure 10. PDD and PDSST data were calculated from daily means derived from the ERA-I reanalysis data. Figure 10. Ratio of inshore to offshore wave power along five distinct sections of Barter Island as a function of offshore incident wave directions. WB, EB, and CB denote the western, eastern, and central bluffs, respectively. Eastern bluff sections are further divided by moderate to very low bluff recession rates (see Figure 7a). Ratios were computed from wave model runs of 2000 individual events of varying wind and offshore wave conditions.
Median wave power values were applied to a continuous hourly hindcast (1979 through 2019) of offshore ERA5 wave conditions to estimate annual cumulative wave power fronting each of the five bluffed coastal sections ( Figure 11). The presence of landfast and pack ice during most of the year limited the analysis period to late May through early December. Annual variations in timing of spring season sea ice breakup in conjunction with wind energy of sufficient magnitude to generate waves can be inferred from the Figure 10. Ratio of inshore to offshore wave power along five distinct sections of Barter Island as a function of offshore incident wave directions. WB, EB, and CB denote the western, eastern, and central bluffs, respectively. Eastern bluff sections are further divided by moderate to very low bluff recession rates (see Figure 7a). Ratios were computed from wave model runs of 2000 individual events of varying wind and offshore wave conditions. flattening, indicating no further wave energy for that particular year. Superimposed upon the individual annual years are solid-colored lines showing the long-term mean, 1979 through 2009, and averages over the latter 10 years split into two 5 year time periods, 2010 through 2014 and 2015 through 2019. The results indicate that inshore wave power has increased four-fold over the past decade compared to the long-term mean, and that the increase is continuing (purple line > blue line), commensurate with the increasing bluff retreat rates.  To evaluate the potential for ambient air and seawater temperatures as drivers of thermal denudation and thawing at the base of the bluff during contact, cumulative positive air and sea surface temperature degree days (PDD and PDSST) were computed ( Figure 11F,G). Because SSTs are indicative of sea ice presence, the spread of the PDSSTs along the monthly x-axis and time at which the curves flatten (indicating cooling and refreeze) were, as expected, similar to the WP plots. The percent change was, however, half that of WP, with~230% increase for the most recent period relative to the long-term mean. Finally, considering the PDD as a potential indicator for thermal denudation, the percent increase relative to the long-term mean was in the order of 160%.
Bluff retreat rates are plotted against cumulative WP, PDD, and PDSST in Figure 12. The physical driver variables (WP, PDD, PDSST) were normalized by the long-term means multiplied by the bluff heights (D i ) to obtain units of meters to match the bluff retreat dimension. Measured bluff heights at each transect and for each time point beginning in 2009 were used when available; bluff heights were assumed to be constant for the preceding years when no data were available. Statistically significant (p-value < 0.005) least-squares linear fits were found for WP (r 2 = 0.92) and PDSST r 2 = 0.90 and only after removal of at least five outlier points for PDD (r 2 = 0.37). (Figure 12A,B). To better characterize the relative contributions of each driver to the measured retreat, a multilinear regression with normalized WP, PDD, and PDSST each multiplied by D i was performed. The resulting slopes were found to be 2.93, 0.17, and 0.04 for WP, PDSST, and PDD, respectively (r 2 = 0.92, Figure 12D). The magnitudes of the slopes indicated that recession was largely driven by mechanical thermoabrasion followed by thermodenudation. Whereas these relationships showed a pattern in an overall sense, the relative contributions of each WP and PDSST to that of PDD was expected to vary significantly across the study area due to nonhomogenous stratigraphy, massive ice content, erratic presence of ice wedges, and evolving surficial gullies. WP, with ~230% increase for the most recent period relative to the long-term mean. F nally, considering the PDD as a potential indicator for thermal denudation, the percen increase relative to the long-term mean was in the order of 160%.
Bluff retreat rates are plotted against cumulative WP, PDD, and PDSST in Figure 1 The physical driver variables (WP, PDD, PDSST) were normalized by the long-term means multiplied by the bluff heights ( ) to obtain units of meters to match the blu retreat dimension. Measured bluff heights at each transect and for each time point begin ning in 2009 were used when available; bluff heights were assumed to be constant for th preceding years when no data were available. Statistically significant (p-value < 0.00 least-squares linear fits were found for WP ( = 0.92) and PDSST ( = 0.90) and on after removal of at least five outlier points for PDD ( = 0.37). (Figure 12A,B). To bette characterize the relative contributions of each driver to the measured retreat, a multilinea regression with normalized WP, PDD, and PDSST each multiplied by was performed The resulting slopes were found to be 2.93, 0.17, and 0.04 for WP, PDSST, and PDD, r spectively ( = 0.92, Figure 12D). The magnitudes of the slopes indicated that recessio was largely driven by mechanical thermoabrasion followed by thermodenudatio Whereas these relationships showed a pattern in an overall sense, the relative contribu tions of each WP and PDSST to that of PDD was expected to vary significantly across th study area due to nonhomogenous stratigraphy, massive ice content, erratic presence o ice wedges, and evolving surficial gullies.

Discussion
The beaches, spits, and barrier islands around Barter Island have continuously retreated to the south and extended to the west through time. Their southerly migration is likely due to a combination of littoral transport, overwash of the islands, and deposition of material into the adjacent lagoons during periods of high waves and elevated water levels. Increasing storm frequency and intensity, combined with sea-level rise, could yield more frequent barrier breaching (breakthroughs) resulting in deterioration of the barrier and reduction in the sheltering capacity provided to the shores of Barter Island. The apparent acceleration in the rate of migration of the western termini of Barter Island and Bernard Spit, beginning sometime between 1987 and 2000 ( Figure 4A), may be related to a change in incident wave directions over a similar time period. This also coincides with the 1987 to mid-to late 2000s period between the predominantly negative (prior to 1987) and positive (mid-2000s to 2020) phases of WP, OWS, and temperature anomalies ( Figure 9). Mean monthly incident wave directions analyzed from the ERA5 offshore wave data (~30 km offshore in~20 m water depth) were primarily from the northeast but appear to have undergone a regime shift, becoming more easterly, between about 1994 and 1999 ( Figure 4B). This easterly shift in incident wave energy likely increased the magnitude of littoral transport toward the west. This potential increase in westward transport and extension of Bernard Spit may also be related to the progradation and widening of the beach fronting the community of Kaktovik and the eastern portion of the Barter Island bluffs, as similarly suggested in Kupilik et al. 2020 [47]. It remains unknown whether the directional change is related to underlying climate patterns such as the Arctic Oscillation or persistent climatic changes. The coastal permafrost bluffs at Barter Island showed relatively modest decadal and interannual erosion rates prior to 2008 (mean: 0.7 to 1.9 m/yr), which were similar to the long-term (~60 year) mean erosional shoreline change rate measured along the Alaskan Beaufort Sea coast (1.8 m/yr; [9]). Beginning in 2008, interannual rates begin to vary widely, from some of the lowest (0.9 m/yr) to some of the highest (6.6 m/yr) mean rates measured during the study period. After 2015, all years exceeded the long-term rate by more than nearly 50% (2.9 to 4.1 m/yr).
Does time-averaging over multiple years in the more limited observational data prior to 2000 mask the natural variability of the highly dynamic coastal environments and dilute the earlier calculated rate of change? Do these observations reflect expected decadal scale variability within the natural system? Alternately, do these data suggest the permafrost bluffs and beaches are responding to larger global or Arctic Ocean scale changes in the forcing physical drivers resulting in a 'new normal'? The answers are probably a combination of these and other factors. Gibbs et al. (2019) [29] investigated 3-dimensional patterns of annual change on the Barter Islands bluffs over a single year that represented typical to slightly below average conditions (2014)(2015). Based on the observed morphologic change and temperature and wave-impact records, they proposed seasonal patterns changed where thermal processes (thermodenudation), primarily indicated by bluff dewatering, surface sloughing, and thaw slumping, and formation of overhanging visors, likely occur throughout the ice-free season but with a relatively larger influence earlier in the season when air temperatures are warmer, daylight is longer, storms are less frequent, and waves rarely reach and erode material from the base of the bluffs. As the ice-free season continues through the summer and into early fall, mechanical and thermo-erosional (thermoabrasional) processes begin to dominate as daylight wanes, temperatures decrease, and numerous large and powerful storms driven by low-pressure storms and easterly winds bring storm surge and waves that erode material from the base of the bluffs and undercut the base through thermo-erosional niching. The loss of the bluff toe and niche development leads to destabilization of the overlying bluff material, which can lead to large block failures, rotational slumps, and surface sloughing ( Figure 13). Debris aprons develop quickly once waves recede from the base of the bluffs, but likely little change occurs to the bluffs or debris aprons once sea-ice formation is initiated and is present to protect the coast. These patterns and high rates of change are believed to be broadly representative of coastal permafrost bluffs found along many high-latitude coastlines worldwide.
Remote Sens. 2021, 13, x FOR PEER REVIEW 21 of 25 erosional niching. The loss of the bluff toe and niche development leads to destabilization of the overlying bluff material, which can lead to large block failures, rotational slumps, and surface sloughing ( Figure 13). Debris aprons develop quickly once waves recede from the base of the bluffs, but likely little change occurs to the bluffs or debris aprons once seaice formation is initiated and is present to protect the coast. These patterns and high rates of change are believed to be broadly representative of coastal permafrost bluffs found along many high-latitude coastlines worldwide. A detailed investigation of the physical drivers of coastal change and morphological response trends are complex and beyond the scope of this paper. However, we began to investigate the possible influence of wave power, direction, air temperature, and water temperature on the rates and patterns of Barter Island bluff retreat.
We hypothesize that the apparent acceleration in retreat rates at Barter Island may be related to increases in both thermodenudation and the number of niche-forming and block-collapsing episodes associated with higher temperatures and longer ice-free conditions in the Beaufort Sea, similar to observations at Drew Point [6].
Comparison between anomalous cumulative annual wave power and temperatures with coastal change rates (Figures 4-6 and 10; Table 2) illustrate the strong influence of these drivers on bluff retreat. The identified long-term increase in wave power and temperatures, particularly after 2003 and 2010, respectively, that also corresponded to a period of high erosion rates at Barter Island, provides a context on decadal-scale climate conditions which can drive coastal change. Further investigation into how broadly these anomalies are recognized along the Arctic coast, as well as potential linkages between the observed anomalies and rates of coastal change with climate indices such as the El Niño Southern Oscillation, Pacific Decadal Oscillation, and/or Arctic Oscillation will help improve the understanding of the interplay between these physical drivers and their influence on coastal change.
Elevated water levels associated with storm surge, high wave energy, and higher air and water temperatures, which are increasing, are all factors that contribute to increased bluff retreat. Teasing apart the relative contribution or importance to annual change from thermo-erosional versus thermomechanical processes is difficult and requires higher frequency driver and response comparisons than included in this paper. However, the volume-loss response due to the two processes are quite different. Thermal processes such as surface sloughing, debris flows, and thaw slumping are typically slow with relatively low volume loss but are temporally consistent. In contrast, thermomechanical processes such A detailed investigation of the physical drivers of coastal change and morphological response trends are complex and beyond the scope of this paper. However, we began to investigate the possible influence of wave power, direction, air temperature, and water temperature on the rates and patterns of Barter Island bluff retreat.
We hypothesize that the apparent acceleration in retreat rates at Barter Island may be related to increases in both thermodenudation and the number of niche-forming and blockcollapsing episodes associated with higher temperatures and longer ice-free conditions in the Beaufort Sea, similar to observations at Drew Point [6].
Comparison between anomalous cumulative annual wave power and temperatures with coastal change rates (Figures 4-6 and 10; Table 2) illustrate the strong influence of these drivers on bluff retreat. The identified long-term increase in wave power and temperatures, particularly after 2003 and 2010, respectively, that also corresponded to a period of high erosion rates at Barter Island, provides a context on decadal-scale climate conditions which can drive coastal change. Further investigation into how broadly these anomalies are recognized along the Arctic coast, as well as potential linkages between the observed anomalies and rates of coastal change with climate indices such as the El Niño Southern Oscillation, Pacific Decadal Oscillation, and/or Arctic Oscillation will help improve the understanding of the interplay between these physical drivers and their influence on coastal change.
Elevated water levels associated with storm surge, high wave energy, and higher air and water temperatures, which are increasing, are all factors that contribute to increased bluff retreat. Teasing apart the relative contribution or importance to annual change from thermo-erosional versus thermomechanical processes is difficult and requires higher frequency driver and response comparisons than included in this paper. However, the volume-loss response due to the two processes are quite different. Thermal processes such as surface sloughing, debris flows, and thaw slumping are typically slow with relatively low volume loss but are temporally consistent. In contrast, thermomechanical processes such as niche formation and block failure are rapid and have the potential to fail very large volumes of material but are episodic. Figure 6 provides a good example of the spatial extent over which block failures can occur (nearly all failure in the central bluff portion of the study area) and the potential volume of material associated with the failures.
The data and methods employed in this study indicate that the recession was dominated by thermomechanical processes, similar to [47]. In that study, the thermomechanical component was related to wave heights of the 8th power, making that component strongly related to wave conditions and less to denudation. In this study, we found that denudation contributed about a quarter of the amount that thermally driven erosion by water contact contributed and <2% compared to mechanical wave erosion.
Only limited data and anecdotal evidence on the relationship between coastal change, storm surge, high wave events, and flooding potential are available for this area. As part of a study to evaluate the stability of Arey Island with respect to sea-level rise and changing storm conditions associated with 21st century climate change, Erikson et al. (2020) [4] derived hindcast (1981 to 2010) and future (2011 to 2100) conditions of the open-water season, wave heights and directions, storm surge, and total water levels using a suite of numerical models. Similar to this study, they noted an expansion in the duration of the open-water season, with nearly equal rates of expansion (1.2 days per year) in the spring and early summer compared to the fall. Their results indicated that the onset and end of the open-water season will shift earlier and later by 20 days per decade and 40 days per decade, respectively, and that the end of the open-water season will extend well into November by the year 2100. Owing in part to the extended open-water season, the number of northwesterly storms affecting the region is expected to increase from less than five events per year prior to 2010 to as many as 30 events per year by the end of the 21st century. Deep water wave heights are projected to increase with increases in storm duration and wind magnitude and peak wave periods, which have implications for increasing water level runup at the shore, potentially increasing wave and water contact with the base of the bluffs.

Conclusions
This study used a rich observational dataset to quantify multidecadal to annual rates and patterns of coastal change at Barter Island, Alaska. The study is unique with respect to the high temporal resolution of the observational data as well as the detailed comparison of bluff change with wave power and sea surface temperatures as indicators of thermoabrasion and ambient air temperatures as proxies for thermodenudation.
The beaches and spits of Barter Island and Bernard Spit showed a generally consistent pattern of southerly (landward) erosion and westerly migration. Beach accretion and progradation was observed only fronting the community of Kaktovik and portions of the LRRS. Based on the historical shoreline change trends and the observed shift to the north east of incident wave energy, continued westerly migration of Bernard Spit is likely, and will continue to shelter to the adjacent mainland leading to continued westerly beach accretion and progradation which will likely protect the bluffs from direct waves, leading to reduced erosion rates for the eastern portion of the bluffs.
Permafrost bluff retreat was continuous over the study period, but rapid retreat was episodic and the rate of change was variable through space and time. Comparison of bluff change rates on decadal time periods showed a steady increase in retreat rates through time, with a near three-fold increase in mean (and maximum) rates of retreat over the last two decades compared to the previous five decades. On annual time scales, episodic erosion driven by thermomechanical undercutting of the bluff base results in some of the highest retreat rates measured along the Beaufort Sea coast of Alaska [6,8,9,29].
Interannual bluff retreat rates showed large temporal and spatial variations, which is characteristic of most coastal environments; however, an overall increase in bluff retreat rates since 2000, and consistently high rates since 2015, suggests acceleration in coastal erosion that is independent of the large spatial and temporal variations observed on an interannual basis.
Similar to the western termini of the Bernard Spit and Barter Island spit, the apex or node of Barter Island, which is generally the same region of highest retreat rates, has also migrated to the west, approximately 800 m since 1950, resulting in an overall straightening of the shoreline, possibly caused by a divergence in longshore transport due to changes in incident wave directions and/or the westerly migration of Bernard Spit.
Rates and patterns of bluff retreat were correlated to incident wave energy and air and water temperatures. Wave energy was found to be the dominant driver, followed by sea surface temperatures, and warming air temperatures and can be considered proxies to evaluate thermo-erosion and denudation. Normalized anomalies of cumulative wave energy, duration of open water, and air and sea temperature show at least three distinct phases: a negative phase prior to 1987, a mixed phase between 1987 and the early to late 2000s, followed by a positive phase extending to 2020. The duration of the open water season has more than tripled since 1979, increasing from approximately 40 to 140 days. Acceleration in retreat rates at Barter Island may be related to increases in both thermodenudation, associated with increasing air temperature, and the number of niche-forming and block-collapsing episodes associated with higher air and water temperature, more frequent storms, and longer ice-free conditions in the Beaufort Sea.
Although this study was focused on a relatively small segment of the Arctic coast, the observed rates, patterns, and drivers of permafrost bluff and barrier change are thought to be broadly representative of similar environments found along high-latitude coastlines worldwide. Continued monitoring along with a more comprehensive evaluation of the processes driving coastal change will improve the overall understanding of these complex coastal systems and their future response to a changing climate.

Data Availability Statement:
The datasets analyzed in this study are available from Gibbs et al., (2020) [40].