Mapping Permafrost Variability and Degradation Using Seismic Surface Waves, Electrical Resistivity, and Temperature Sensing: A Case Study in Arctic Alaska

Subsurface processes significantly influence surface dynamics in permafrost regions, necessitating utilizing diverse geophysical methods to reliably constrain permafrost characteristics. This research uses multiple geophysical techniques to explore the spatial variability of permafrost in undisturbed tundra and its degradation in disturbed tundra in Utqiaġvik, Alaska. Here, we integrate multiple quantitative techniques, including multichannel analysis of surface waves (MASW), electrical resistivity tomography (ERT), and ground temperature sensing, to study heterogeneity in permafrost’s geophysical characteristics. MASW results reveal active layer shear wave velocities (Vs) between 240 and 370 m/s, and permafrost Vs between 450 and 1,700 m/s, typically showing a low‐high‐low velocity pattern. Additionally, we find an inverse relationship between in situ Vs and ground temperature measurements. The Vs profiles along with electrical resistivity profiles reveal cryostructures such as cryopeg and ice‐rich zones in the permafrost layer. The integrated results of MASW and ERT provide valuable information for characterizing permafrost heterogeneity and cryostructure. Corroboration of these geophysical observations with permafrost core samples’ stratigraphies and salinity measurements further validates these findings. This combination of geophysical and temperature sensing methods along with permafrost core sampling confirms a robust approach for assessing permafrost’s spatial variability in coastal environments. Our results also indicate that civil infrastructure systems such as gravel roads and pile foundations affect permafrost by thickening the active layer, lowering the Vs, and reducing heterogeneity. We show how the resulting Vs profiles can be used to estimate key parameters for designing buildings in permafrost regions and maintaining existing infrastructure in polar regions.


Global Warming Impacts on Permafrost Degradation and Ground Thermal Regime
Permafrost refers to any ground material that remains below 0°C for at least two consecutive years.It is distributed across approximately 25% of the land surface in the northern hemisphere and is highly sensitive to atmospheric temperature variation primarily caused by global warming (Biskaborn et al., 2019;Lantuit et al., 2012).The Arctic average annual surface temperature has increased by 3.1°C from 1971 to 2019, which is three times faster than the global rate (AMAP, 2021;IPCC, 2021;Rantanen et al., 2022).In Utqiaġvik (formerly known as Barrow), Alaska (study area), the average annual air temperature has risen over 4°C since 1980, and Arctic Alaska's permafrost has warmed by 1-3°C (Nicolsky et al., 2017;Thoman & Walsh, 2019).Increasing temperature in the high-latitude permafrost regions leads to permafrost degradation, which includes permafrost warming, active layer thickening, and thaw-related hazards such as the development of taliks, ground subsidence and thermokarst in low-lying areas, mass wasting on slopes, and thermal erosion and abrasion along riverbanks and coasts (Hjort et al., 2022).Global warming causes contaminated and industrial sites in regions of stable permafrost to thaw, posing a significant environmental threat (Langer et al., 2023).Permafrost degradation drives serious changes in local geomorphology, hydrology, vegetation, wildlife dynamics, and greenhouse gas emissions (Hjort et al., 2022;Streletskiy et al., 2015).
Permafrost research often focuses primarily on ground temperature due to its direct effect on physical and biogeochemical soil processes, but permafrost is also affected by air temperature, snow cover, soil moisture, vegetation cover, and soil properties (Lantuit et al., 2012;Smith et al., 2022).In situ monitoring using thermistors and thermocouples has shown that permafrost temperatures are increasing, leading to thawing and degradation (Biskaborn et al., 2019;Nicolsky et al., 2009Nicolsky et al., , 2017;;Romanovsky et al., 2010;Shiklomanov et al., 2010).Understanding the ground's thermal state in permafrost regions is crucial to model and mitigate climate change impacts.
Understanding the spatial heterogeneity of permafrost in Arctic tundra is important for studying geomorphological and ecosystem variations under climate change as well as potential engineering impacts.Permafrost structure is often complex due to fine-scale spatial heterogeneity of properties such as temperature and ice content.Temperature, water saturation, and ice content influence seismic wave velocities, including shear wave velocity (V s ) and compressional wave velocity (V P ) (Coduto, 1999;Hjort et al., 2022;Ji et al., 2023;Liew et al., 2022;Rocha dos Santos et al., 2022).Salinity is a key factor in coastal permafrost environments such as Utqiagvik, Alaska.Studies by Brown (1969) and O'Sullivan (1966) reveal how salinity affects permafrost's geochemistry and physical properties, such as its freezing point and stability.Dafflon et al.'s research (2016Dafflon et al.'s research ( , 2017) ) further demonstrates the impact of salinity on the distribution and characteristics of shallow permafrost, highlighting its connection with vegetation patterns and soil properties.Meyer et al. (2010) provide the historical context, showing how historical salinity variations have influenced the permafrost landscape over time.
Here, we provide a brief overview of the active layer and some of the key permafrost structures investigated in this paper, including cryopegs, ice-rich zones (e.g., lenses), and thermokarst lakes.Jafarov et al. (2018) showed that the active layer thickness (ALT) of undisturbed tundra (without infrastructure development and low to medium water content) near Elson Lagoon in Utqiaġvik, Alaska is approximately 0.2-0.6 m, measured in August 2013.A thermokarst is formed when the thermal equilibrium shifts, allowing the ground ice to thaw.Talik is a layer or body of year-round unfrozen ground (usually above 0°C) occurring in a permafrost zone due to a local anomaly in thermal, hydrological, or hydrochemical conditions (e.g., underneath thermokarst lakes and rivers).Cryopegs can have temperatures below 0°C, but freezing is prevented by freezing-point depression due to the dissolved-solid content of the pore water (Van Everdingen, 1998).As the climate warms, the annual ground temperature increases, and annual thawing deepens until a certain threshold is met, after which a talik develops.Ice wedges form when water seeps into cracks in the ground during summer and then freezes during winter.The distribution of ice formations and ice content within the permafrost layer is highly variable (Liu et al., 2021).Previous studies confirm that the ice content in the permafrost around Barrow is very high in the upper part and decreases with depth (Brown et al., 1980).For the coastal plain along the Beaufort Sea from Point Barrow to the Canadian border, Kanevskiy et al., 2013 reported an average total ice content (i.e., ice wedge, segregated, and pore ice) of 83% and 82% for the primary surface and the drained-lake basins, respectively.This results in a landscape vulnerable to widespread subsidence and thermokarst development, the magnitude of which may vary widely depending on surficial geology, ground ice volume, and the extent of past thermokarst activity (Farquharson et al., 2019).The warming and thawing of ice-rich permafrost pose changes in its interactions within the built environment (Hjort et al., 2022).
Permafrost exhibits vastly variable properties between thawed and frozen states due to the phase change of water, impacting its strength and bearing capacity, which can lead to infrastructure failure (Hjort et al., 2022).The interaction between permafrost and civil infrastructure contributes to permafrost degradation and increases construction and maintenance costs (Streletskiy et al., 2012(Streletskiy et al., , 2015)).While existing research primarily focuses on the influence of degrading permafrost on infrastructure, it is crucial to consider the impact of civil infrastructure on permafrost as well.Different foundations and architecture, such as pile foundations and gravel roads, introduce thermal and physical impacts that can disturb the natural environment and alter adjacent tundra ecosystems (Walker et al., 2022).As climate change continues, the vulnerability of both civil infrastructure and permafrost systems increases, necessitating detailed knowledge of risk exposure in current and future infrastructure areas (Hjort et al., 2022;Melvin et al., 2017).Understanding the influence of civil infrastructure on degrading permafrost allows for a realistic risk assessment.

Seismic and Geoelectrical Methods in Permafrost Regions
Seismic imaging is a commonly used technique for characterizing the subsurface in permafrost regions (e.g., Justice & Zuba, 1986;Miller et al., 2000;Ramachandran et al., 2011), because seismic wave velocities, including shear wave velocity and compressional wave velocity, are sensitive to temperature, water saturation, and ice content (Coduto, 1999;Hjort et al., 2022;Ji et al., 2023;Liew et al., 2022).One of the main advantages of seismic imaging in permafrost regions is its ability to provide detailed information about the distribution and continuity of permafrost and the nature of the underlying soils, even 3D profiles with high resolution (e.g., Ramachandran et al., 2011;Schwamborn et al., 2002).This information is important for various applications, such as infrastructure planning and design, resource exploration, and environmental monitoring.Seismic refraction is a surface geophysics method that utilizes the refraction of body waves through layered media (Scott &

Study Area
Permafrost zones underlie 80% of Alaska, including 29% continuous permafrost (Jorgenson et al., 2008).The North Slope Borough is entirely within the continuous permafrost zone (Ferrians, 1965;Kerkering, 2008), as shown in Figure 1a.The permafrost in Utqiaġvik, Alaska, is continuous and has a thickness of approximately 200-400 m (Jorgenson et al., 2008).Elson Lagoon forms the eastern land boundary of the study area, as shown in Figure 1b.The ALT of undisturbed tundra near Elson Lagoon is approximately 0.2-0.6 m, and the soil volumetric water content varies from 17% to 88%, measured in August 2013 (Jafarov et al., 2018).The average ALT of the study area on the tundra in Utqiaġvik, Alaska, is also shallow (roughly less than 1.0 m) consisting of three distinct layers: the acrotelm (top), the catotelm (middle), and the mineral soil (bottom) (Chen et al., 2020).The ground conditions vary from dry to marshy, with surface vegetation.The seismic survey (MASW) operations were performed on 6-12 August 2022.The plan view showing the layout of seismic survey lines and the temperature measurement locations is presented in Figure 1.

The MASW Technique
Surface waves can be generated by an active source, such as a hammer, weight drop, or vibroseis, or by a passively recorded source, such as anthropogenic, traffic, or several other environmental sources (e.g., ocean waves, wind), and these waves are recorded by an array of geophones.Because surface waves can provide information on the subsurface velocities over a wide range of frequencies and wavelengths, the MASW technique can generate high-quality velocity models.MASW is often used to produce 1D velocity profiles, and it can also be used to assess the lateral variability of the subsurface shear wave velocities, which is essential for characterizing subsurface heterogeneity and identifying areas of potential geotechnical concern.The MASW method can be applied in a wide range of geological environments and has the advantages of being non-invasive, cost-effective, and capable of providing high-resolution V s profiles to depths of up to several tens of meters.
In active-source surface seismic surveys, over two-thirds of the total seismic energy generated by compressional waves is transmitted to Rayleigh waves, sometimes referred to as "ground roll" (Park et al., 1999b).An example of surface waves in our collected data is shown in Figure 2a.Surface wave energy decays exponentially with depth beneath the surface.Longer wavelength (i.e., longer-period and lower-frequency) surface waves travel deeper, thus containing more information about deeper parts of the study area.Shorter wavelength (i.e., shorterperiod and higher-frequency) surface waves travel shallower, thus containing more information about shallower parts of the study area.Surface waves are dispersive, meaning each wavelength propagates at different phase velocities in a layered medium.Thus, we can analyze the phase velocity of different frequency bands (corresponding to different wavelengths) and estimate the velocity profile of the subsurface.
Rayleigh wave dispersion curves describe the velocity at which each wavelength travels.To determine Rayleigh wave dispersion curves, we use the phase shift method, which provides accurate fundamental-mode phase velocities even when only four geophones are used (Dal Moro et al., 2003;Park et al., 1999b).Our detailed procedure is provided in Supporting Information S1. Figure 2 represents an example of the inversion procedure for estimating a 1D velocity profile.First, we applied a 7.50-327.68Hz bandpass filter to all traces in a shot gather to remove high-frequency noise, and we muted noisy traces (Figure 2a).Then, dispersion images are calculated to determine phase velocities at each frequency (Figure 2b).The fundamental mode velocity profile of the Rayleigh surface wave are extracted from the dispersion image (Figure 2c).After wavelength-depth conversion, we generated an initial model based on the phase velocity picks.Finally, a non-linear least squares inversion method (Xia et al., 1999) was applied to the dispersion curve to reconstruct the V s velocity model (Figure 2d) using the SeisImagerSW software.The minimal depth at which shear wave velocity can be reliably inferred through inversion is contingent on a confluence of site-specific variables.This interpretable depth is not a fixed measure but is instead modulated by an array of factors, including the fidelity of seismic data, particularly its high-frequency components, and the proximal spacing relative to the seismic source.Acceptable final 1D models should have a root mean square (RMS) error of the difference between the theoretical dispersion and measured dispersion curves (Figure 2c) below 5% (SeisI-magerSW TM Manual, 2009).With this procedure, we obtain an average 1D V s vertical profile along the seismic line.
For 2D V s estimation, we carry out the same pre-processing and then perform dispersion analysis using the common mid-point (CMP) cross-correlation gathers.The CMP cross-correlation method can increase the signalto-noise ratio of the dispersion spectrum and improve lateral resolution (Hayashi & Suzuki, 2004).Here, we briefly discuss the following steps for CMP cross-correlation analysis: First, the cross-correlation functions between every pair of traces in each shot gather are calculated.Second, correlation traces with a common midpoint are collected, and those with the same spacing are stacked.The resultant cross-correlation gathers resemble shot gathers and are known as CMP cross-correlation gathers.Third, dispersion image of surface waves are calculated using the MASW technique (Park et al., 1999b) from the CMP cross-correlation gathers.Finally, a non-linear least squares inversion method is applied (Xia et al., 1999) to the dispersion curves for reconstructing a 2D V s model.As a general guideline, acceptable 2D models should result in an RMS below 15% (SeisIma-gerSW TM Manual, 2009).
MASW is limited to shallow depth investigations, typically up to 30 m.Beyond this depth, the resolution and accuracy of the method decrease significantly.Moreover, MASW is more suitable for homogeneous soil conditions and struggles to characterize laterally varying structures or complex geological settings accurately (Boiero & Socco, 2011;Evangelista & Santucci de Magistris, 2015).To overcome these limitations in the future, seismic refraction or reflection could be employed over a larger region as a complementary technique, although reliable picking of the un-aliased arrivals in this complex near-surface zone will likely require more than 24 geophones and a high-frequency seismic source.The depth and compression velocity information of different subsurface layers can be determined by analyzing the travel time data.Seismic refraction has the advantage of investigating deeper depths, making it useful for studying subsurface structures beyond the reach of MASW.The material in Text S2 of the Supporting Information S1 can serve as a basis for a future permafrost study combining MASW and seismic refraction.Advances in data acquisition and processing should be taken into account to ensure the best possible outcomes in future investigations.By integrating the strengths of these two techniques, it may be possible to enhance the accuracy and depth range of permafrost characterization.

Electrical Resistivity Tomography Method
Surface electrical resistivity surveying is based on the principle that the distribution of electrical potential in the ground around a current-carrying electrode depends on the electrical resistivities and distribution of the surrounding soils and rocks.The usual practice in the field is to apply an electrical direct current (DC) or alternating current (AC) of low frequency.The voltage between two potential electrodes and the current between two other electrodes are measured during the ERT survey.Measurements were provided using various pairs of electrodes along the transect.Increasing the spacing between electrodes allows for deeper investigation depth.Therefore, apparent resistivity values can be obtained through the voltage, current, and geometry of array electrodes for each point laterally and with depth.The measured apparent resistivity values are used in the inversion to obtain resistivity models for subsurface materials.
After data acquisition, the data was checked for large measurement errors (>2%) using Prosys II software (IRIS Instruments), which was provided with the Syscal equipment.During the data acquisition process, a roll-along technique was implemented, involving the relocation of the array's center to a distance of 18 m (approximately half of the initial length equal to 35.5 m).To assess the quality of the data, we computed the errors between measurements taken after this center movement.Approximately 20% of the total measurements were repeated, and the average error was found to be less than 1%.These metrics provide valuable insights into the overall data quality, indicating a relatively low average error and the redundancy of about 20% in the measurements.For the purpose of joint inversion, we combined Wenner-Schlumberger (WS) and Dipole-Dipole (DD) data into one data set.
Res2Dinv software (Res2Dinv Manual, 2006) was used to generate a resistivity model from apparent resistivity values.The inversion software uses a smoothness-constrained least squares method with L2-norm (Loke & Barker, 1996).For the model mesh, we used half-electrode spacing (0.25 m) in the horizontal direction.In the vertical direction, the thickness for the first layer was set to 0.125 m with a 1.1 factor to increase thickness with depth.
The inverted resistivity model can also be influenced by topographic changes along a transect, particularly on the top cells.Minor variations in topography may have negligible effects on the inversion process, but substantial changes can introduce distortions.Elevation changes along the transect are generally insignificant and less than 0.5 m on the edges of ice wedge polygons for our study site (at MASW 3 location).As a result, topography has not been integrated into the inversion process.In some cases, degraded polygons have large elevation changes, and it is crucial to carefully evaluate the significance of topographic changes for these areas when dealing with geophysical inversion models.
We used a robust constraint and defined an initial half-space resistivity value of 10 Ωm to trace spatial cryopegs distribution based on low resistivities observed in previous studies (Hubbard et al., 2013;Overduin et al., 2012;Yoshikawa et al., 2004).A 2D ERT model was generated to invert the WS and DD arrays jointly, and the result achieved an RMS error of 1.8%.Across different iterations, specifically the 3rd, 4th, and 5th, the percentage changes were 5.8%, 2.8%, and 1.86%, respectively.When comparing the results of inverted resistivity with varying RMS errors, notable variations were observed primarily in the lower cells of the model.The overall changes between the 4th and 5th iterations were evident in both resistivity and the bottom boundary for the deepest contrast layer.Nonetheless, these changes remained generally stable and did not substantially affect our interpretation results.
ERT and other electromagnetic methods provide inherently non-unique solutions for mapped reconstructions of subsurface electrical properties.This non-uniqueness means that the measured data can be explained equally well by multiple models.To improve the reliability of subsurface interpretations, we use soil core sampling and comparison with other geophysical methods, such as the MASW technique.

Correlation of Soil Mechanical Properties With Shear Wave Velocity
Soil's mechanical properties can be determined using V s , which is a commonly used geotechnical and geophysical parameter.There are empirical or analytical correlations between V s and several other soil properties.For instance, there is a positive correlation between V s and soil stiffness parameters such as shear modulus (G) and elastic modulus (Young modulus, E).Stiffer soils generally exhibit higher shear wave velocities and are correlated with higher soil strength parameters such as undrained shear strength (S u ) and peak shear strength.In addition, V s is inversely correlated with soil porosity, where lower V s values are often observed in soils with higher porosity.For a given soil type, with constant parameters such as temperature and ice content, Vs is directly proportional to soil density, indicating that denser soils generally exhibit higher V s values.Furthermore, soil classification, which determines the foundation design, directly correlates with V S30 and S u30 , which are the averages of V s and S u in the top 30 m (ASCE/SEI 7-16, 2017).Therefore, understanding the soil classification and V s is crucial for assessing the seismic performance and stability of foundations and for designing appropriate foundation systems.
For soil that is elastic, isotropic, and homogeneous, the elastic theory can be used to establish the following relationship between elastic modulus and seismic wave velocity: (1) ,and (2) where, μ is the Poisson's ratio, γ is the unit weight of the media, and g is the gravitational acceleration, which can affect the foundation design in various ways such as foundation type, foundation settlement, and allowable vertical and lateral loads (Coduto, 1999).The relationship between elastic modulus and bearing capacity of soils could depend on several factors such as soil classification, stress history, foundation type, etc.The elastic settlement beneath a flexible footing placed on the ground surface can be calculated as follows (Terzaghi et al., 1996): where q is the surcharge load, B is the width of the footing, and I f is the influence factor, which is a function of the ratio of length to width ( L B ) and the thickness of the compressible layer.We will highlight the importance of V s monitoring on changes in elastic modulus and elastic settlement in Section 5.6.

Data Acquisition
There are eight seismic survey locations using MASW, six temperature measurement locations using thermistors, one ERT survey location (at MASW 3 location), and five core sampling locations (at Roadside and MASW 1-4 locations), as shown in Figure 1.The coordinates of these locations are provided in Table S1 of the Supporting Information S1.The seismic surveys cover various soil conditions, including disturbed and undisturbed areas.Four seismic surveys (MASW 1-4) were performed on undisturbed tundra permafrost and four (Roadside and NOAA 1-3) on disturbed permafrost, shown in Figure 1c.The seismic surveys performed on disturbed permafrost include one survey along the gravel road near the National Oceanic and Atmospheric Administration (NOAA) facility (Roadside), one survey under the NOAA building (NOAA 1), one survey on the pre-existing building foundation next to the NOAA building (NOAA 2), and one survey on the tundra near the pile foundations (NOAA 3).At the NOAA 2 location, a building was demolished and removed 1 year prior to the seismic survey, but the pile foundations remain in the ground.The core sampling was performed at approximately 1 m from the seismic survey locations using a hand-held sampling drill (boring auger).The sampling depth is up to 103 cm.

Surface Wave Data Acquisition
Each seismic line consists of 24 vertical 4.5 Hz geophones (a 24-channel Geometrics Geode seismograph) positioned on the ground surface.Straight-line seismic profiles have geophone spacing equal to 1 m, which gives us 23 m spread in total.We generated seismic signals using a sledgehammer adjacent to the geophones as well as an extra shot at 5 m offset from the beginning of the lines.The seismic record length was 128 ms with a sample interval of 0.25 ms, and each recording was initiated by a trigger attached to the sledgehammer.No preacquisition filter was used on the seismic data.Note that with the vertical source and the vertical receivers, the type of surface waves we acquired are Rayleigh waves.An example of the collected seismic traces is shown in Figure 2a.

Temperature Data Acquisition
Small holes with 0.02-m diameter were punched in the ground to a depth of 1.5 m in August 2021.Four HOBO TMC6-HD temperature sensors were then lowered into the ground using wooden rods to the depths of 0.02, 0.2, 0.5, and 1.5 m below the surface.At each location, temperature sensors were connected to two 2-channel HOBO U23-003 loggers in September 2022 (Soil temperature measurement information is presented in Figure S3 and Table S2 of Supporting Information S1).The operating temperature range for loggers is 40-100°C with an accuracy of 0.4 and 0.2°C below and above 0°C, respectively.The resolution is 0. prohibited in this study area, reducing the vulnerability of the sensor installation to damage by wildlife animals is a challenging task and several temperature sensor cables were severed by Arctic foxes.Nevertheless, temperature records were collected at the following six sites.Site conditions were described during installation as follows.The Temp 1 profiler was placed in a shallow pond with 10 cm of standing water; The Temp 2 profiler is located near a rim of the flat-center ice-wedge polygon; the site conditions could be described as moist.The Temp 3 station is placed at the center of a high-center polygon with a dry ground surface.The Temp 4 site is almost saturated, with a thin layer of standing water in the middle of the low-center polygon.The Temp 5 profiler is at the rim of the lowcenter polygon, with the dry ground around it.Finally, the Temp 6 profiler is approximately 9 m from shore, where the ground is rather moist.

Electrical Resistivity Tomography Data Acquisition
An electrical resistivity tomography (ERT) survey was conducted to provide in situ resistivity measurements (Ωm) along the MASW 3 transect from 13 September 2022.The ERT station Syscal-Pro 72 (IRIS instruments) and steel electrodes were used to acquire data.Electrodes were placed along the transect using a measuring tape with a 0.5 m spacing.Inverse WS and DD arrays were applied for measurements.The minimum/maximum half electrode spacing was 0.75/17.25 m for WS and 0.5/13.5 m for the DD array.A 50 V output voltage and 250 ms pulse duration were applied during the survey.Contact resistances for most of the electrodes were no more than 1 kΩ.Measurements with errors exceeding 2% were removed during processing.

Borehole Core Sample Acquisition
Permafrost core samples were collected using a hand-held boring auger at MASW1-4 and Roadside locations.The Auger boring technique is particularly effective in permafrost areas due to its ability to efficiently penetrate the frozen ground with limited disturbance to the surrounding soil.Core samples with diameters of 3.8 cm and depths of up to 103 cm were collected to determine active layer depth, identify soil type and density within the active layer and permafrost, evaluate salinity levels, and construct stratigraphy plots.Core samples were kept frozen during the transportation to the University of Alaska Fairbanks for identification and density measurements.The thawed samples were then sent to the Pennsylvania State University for salinity measurements.

Results and Discussion
This section provides an overview of the V s , electrical resistivity, and temperature results.We highlight key findings related to permafrost, including site characterization (Section 5.1), shear wave velocity interpretation (Sections 5.2-5.4),civil infrastructure's influence on permafrost (Section 5.5), and applications to engineering properties and infrastructure design (Section 5.6).

Site Characterization of Disturbed and Undisturbed Permafrost
We calculated 1D and 2D V s models for 8 survey lines after performing inversion on the calculated dispersion curves using the MASW method (Park et al., 1999a(Park et al., , 1999b)).This included both disturbed and undisturbed permafrost regions.Figure 3 shows the 1D V s profiles for the survey locations on the tundra, which are categorized as undisturbed permafrost locations.Figure 3a shows the V s profile for the first location of the MASW survey in the tundra (MASW 1), approximately 500 m from the road and NOAA facility.A similar V s profile is observed in Figures 3b-3d (MASW 2-4) except that the relative high-velocity zone is located relatively at a higher depths than MASW 1. Also, MASW 1 has a more consistent velocity profile at higher depths than MASW 2-4.These could be an effect of the different geology, vegetation, or the effect of anthropogenic activities over the years, as we observed many marked points for previous studies and tracks from vehicles in the field.Figure 3b represents the V s for the second location in the tundra (MASW 2), approximately 1 km from the infrastructure.We observe a very high-velocity zone at a depth of 2-8 m below the surface, representing either a stiff lithology layer or an ice-rich zone.The highest V s layers are located at a depth of 5 m and are as high as 1,700 m/s.This location clearly shows the undisturbed permafrost area with higher V s and higher ice content.Figure 3c illustrates the velocity model for the third location in the tundra (MASW 3), roughly 1.5 km from the road and NOAA building.While the low-high-low V s pattern is obvious, the highest velocity is 1,575 m/s, which is lower than that at the MASW 2 at 1 km.In the field, we observed that as we get further to the tundra, the ground gets wetter as indicated by many ponds in the area.This can also be seen on the satellite map in Figures 8c and 8d (presented later in this paper), where the last two lines (MASW 3 and 4) are located in darker areas that represent higher surface water content.Figure 3d illustrates the velocity model for the last location in the tundra (MASW 4) at roughly 2 km from the road and NOAA building.Like all other locations in the tundra, we observed the low-high-low velocity profile, but the highest velocity zone is located at deeper depths of 7 m.In addition, the highest V s is 1,150 m/s at MASW 4, which is lower than those at the two previous locations in undisturbed permafrost zones (MASW 2 and 3).These last two locations (MASW 3 and 4) represent the effects of vegetation and high surface water content.The previously observed low-high-low velocity profile below the surface by Dou and Ajo-Franklin ( 2014) is captured at all tundra locations (i.e., MASW 1-4 and NOAA 3).
Figures 3e-3g show the 1D V s profiles for the survey locations near the NOAA building, including the roadside, within 1.0 m of the NOAA building (NOAA 1), and ∼80 m away from the NOAA building (NOAA 3), where we expect high disturbance in permafrost.These locations are categorized as disturbed permafrost locations.The V s at disturbed permafrost locations is lower than that at undisturbed permafrost locations, as the vegetation and thus albedo are affected by human activities, and the ice content is lower than that in the undisturbed permafrost Journal of Geophysical Research: Earth Surface 10.1029/2023JF007352 TOUREI ET AL. region.Figure 3e represents the V s for the location within 1.0 m of the gravel road leading to the NOAA facility, where we expect high disturbance in permafrost.The 1D velocity profile represents a low-high-low velocity pattern with the highest V s of 850 m/s.It is noted that the data at the pre-existing demolished building foundation next to the NOAA building (NOAA 2) are low-quality, which resulted in a higher RMS error than the acceptable error.As we performed sledgehammer shots on top of a pre-existing building foundation, the contrast in soil and pile material properties and the resulting scattered energy likely generated relatively larger errors than in other locations.
Figure 4 shows the 2D V s models in the undisturbed (Figures 4a-4d) and disturbed (Figures 4e-4g) permafrost regions similar to Figure 3.A low-high-low velocity pattern is evident in all models, indicating the active layer and transitional zone (with low velocity), ice-rich permafrost (with high velocity), and partially frozen permafrost with scattered ice (low velocity).Previous studies have shown a strong correlation between permafrost temperature and V s (Ji et al., 2023;Kurfurst, 1976;Nakano et al., 1972).The low-high-low velocity pattern is consistent with the general trend of permafrost temperature variation with depths, as in a previous study (Smith et al., 2022).Although the velocity and depth vary spatially, the low-high-low pattern is consistent among all locations.The 2D models capture the spatial variability of permafrost, demonstrating the importance of multichannel seismic surveys and 2D modeling.We will discuss these models in detail in Sections 5.2 and 5.3.Meanwhile, we invert the electrical resistivity model (Figure 5a) from the ERT data, which can delineate several zones compared to the V s model from the MASW 3 transect (Figure 5b).Each zone is characterized by different thicknesses, electrical resistivity, and velocity values.Zone A is characterized by relatively low resistivity up to 200-300 Ωm from the surface down to approximately 0.4 m.This layer represents the active layer and is characterized by low V s of 240 m/s.Zone B is characterized by high resistivity values of 400-2,000 Ωm.The thickness of this layer varies from 1 m at distances of 6-24-4 m at distances 2-6 m.According to drilling data, the zone between 1.4 and 1.7 m is a transition zone (black box marked in the borehole in Figure 5) between a frozen and unfrozen state and represents a boundary between zones B and C. Zone C is characterized by soils with low resistivity of 10-20 Ωm associated with cryopegs development in the study area.The thickness of the layer is about 2.5 m along the profile and decreases at the beginning of the transect.At a depth of 0.5-0.6 m, the salinity measured approximately 0.84 ppt (parts per thousand, roughly equivalent to grams per liter).In contrast, at a deeper level, specifically 1.5-1.6 m, the salinity was significantly higher, measuring 8.03 ppt.Consequently, one can infer that salinity gradually increases from Zone B to Zone C, contributing to the low resistivity observed for cryopegs in Zone C. The characteristics of Zone C, exhibiting low resistivity and relatively high velocity, may be attributed to two potential factors.First, it is plausible that this zone remains frozen, with its temperature lingering below the freezing point of cryopegs.Alternatively, the distinct lithology and soil mineral composition within this zone could be significantly influencing the resistivity and seismic velocity properties of unfrozen sediments in the near-surface, as highlighted in studies by Rossi et al. (2022) and Accaino et al. (2023), respectively.These findings are aligned with previous studies demonstrating the complex and heterogeneous nature of permafrost in the study site (Brown, 1969;Dafflon et al., 2016).Zone D is characterized by relatively high resistivity values of 100-600 Ωm and a thickness of approximately 2 m and is located between two zones of low resistivity (both zones denoted as Zone C).The zone is also characterized by high V s up to 1,600 m/s, typical for frozen material and increased ice content.Due to the relatively high-velocity values, the high-resistivity layer D, and the presence of salt pockets even in ice-rich conditions in Utqiagvik (Iwahana et al., 2021), we interpret this layer as an ice-rich layer.High-velocity values for zone D are also supported by the intra-ice brine pockets, where cryopeg brine is bounded by ice, and is generally in solidified form (Iwahana et al., 2021).
Due to the high contrast of the resistivity of different units, the ERT method helps to identify multiple layers in the upper part of the cross-section that cannot be clearly distinguished using the MASW method due to the lack of high-frequency signals.The MASW results provide useful information to verify ERT results at greater depths.The MASW method reveals a high degree of heterogeneity in the permafrost, possibly due to increased salinity at greater depths, resulting in unfrozen zones.Additionally, integrating ERT and MASW provides a more comprehensive assessment, offering insights into critical aspects, such as the salinity and ice content within the permafrost, which are crucial for understanding its physical properties and behavior.This methodological framework, effective in our specific study area, has potential applicability in diverse permafrost environments, indicating its potential for wider application and enhancing the knowledge of near-surface permafrost dynamics.

Identification of Active Layer Thickness
Significant variations in ALT exist between different landscape types, reflecting the influence of vegetation, substrate, microtopography, and especially soil moisture (Shiklomanov et al., 2010).From Vs profiles shown in Figure 4, we can identify the ALT range (roughly 0.3 m) and the shear wave velocities of the active layer (240-370 m/s) in most locations.However, in some undisturbed permafrost regions, ALT was found to be highly spatially heterogeneous due to differences in subsurface characteristics based on 2D V s profiles shown in Figure 4. Therefore, ALT at some locations may be 0.5-0.6 m, which is consistent with the ALT range estimated from nearby temperature measurements (described in Section 4.2) and the earlier estimation by Jafarov et al. (2018).In contrast, ALT in disturbed permafrost regions presents higher values (0.5-1.0 m) and less spatial heterogeneity.The higher ALT indicates that the ground temperature is slightly higher for the permafrost with human activities than in the undisturbed permafrost.For NOAA 1 (Figure 3f), the high consistency in ALT is because the ground surface is under the NOAA facility, and the topsoil is gravel, which is different from all other locations.The coverage of the building, which produces continuous heat, and the high thermal conductivity of gravel compared to fine-grained soil like peat or silt are likely contributing factors to the temperature consistency of the tested line.The MASW results did not reveal the top shallow active layer of MASW 3 due to the small ALT (0.24 m based on soil sampling) and lack of high-frequency source signal required to image shallow depths.Figure 6 illustrates the stratigraphy of deposits up to 1.0 m at the five sampling locations.The permafrost mainly consists of silty soil with organics with an average unit weight of 19.62 kN/m 3 .Ice layers are evident in the boreholes of MASWs 1 and 2. Based on the boring samples collected at the same sites of the seismic surveys, ALT in the study regions was 0.14-0.28m during sampling.The ALT range identified by MASW in most locations is within a reasonable range compared with the stratigraphy of deposits.

Identification of Spatial Heterogeneity of Permafrost
Spatial heterogeneity within the permafrost layer can be observed and quantified by analyzing V s profiles, including ice-rich permafrost, low-velocity zones, and cryopeg.Shear wave velocities within the permafrost layer range from 450 to 1,700 m/s.The shear wave velocities ice-rich permafrost zones (MASW 2-4) are in the range of 700-1,700 m/s, which is higher than the range of 500-900 m/s from other permafrost locations.Ice-rich permafrost can be identified in 2D V s profiles, such as in Figure 4b (MASW 2), with a high-velocity zone (from 9 to 22 m) with a V s range of 1,300-1,700 m/s.The theoretical V s of pure ice is approximately 1,900 m/s at a temperature near 10°C (Kohnen, 1974).Given that the effective V s of a medium is a weighted average of the components of that material, the regions of the subsurface with velocities near 1,700 m/s are expected to be icerich materials.This indicates that the center area of the ice-rich zone is likely composed primarily of ice layers.However, the gradual increase in velocity near the ice-rich zone at MASW 2 indicates suspended soil around the ice layers.
Figure 1c shows that MASW 3 is located in a wetter area (darker image color is related to higher surface water content), which may lead to open talik regions around the large water body.A potential reason for the talik or cryopeg layer at MASW 3 is salinity, as higher salinity layers exist at MASW 3 due to the proximity of the nearby saline thermokarst lake.This is consistent with our observation of core sample's salinity at MASW 3 (discussed in Section 5.1).

Impacts of Ground Temperature and Ice Structure on Shear Wave Velocities
In Sections 5.4 and 5.5, we focus on the impacts of multiple factors on shear wave velocity in undisturbed permafrost to better understand permafrost behavior and stability.Figure 7 shows the composed profiles of V s in undisturbed permafrost, temperature variation, and cryostratigraphy versus depth.Generally, the rate of temperature decrease lessens with increasing depth.The temperature measurements are derived from several locations in North Slope Borough, Alaska, mainly near the study site.The detailed location and record date of the temperature measurement are presented in Figure S2 and Table S2 of Supporting Information S1.The temperature reveals that the ALT is around 0.5 m, which agrees with the ALT determined from the MASW surveys.
We observed that the depth variation of V s (Figure 7) exhibits a consistent trend across different testing locations near Utqiaġvik.The V s is low (∼250-510 m/s) in the active layer and increases in the permafrost layer to around 1,200-1,700 m/s as depth increases to 5-8 m.Beyond this depth, the shear wave velocities decrease, forming a low-velocity permafrost zone (∼500-700 m/s) beneath the high-velocity permafrost layer.Dou and Ajo-Franklin's (2014) study also reported the existence of a low-velocity permafrost zone at a location approximately 5 km south of the study area.This suggests that low-velocity permafrost zones may exist under the tundra near Elson Lagoon and east of Utqiaġvik.As shown in Figure 7, the V s profiles are correlated with temperature profiles, with higher ground temperature corresponding to lower shear wave velocities of permafrost.We also observe that Vs below 8 m decreases, possibly due to discontinuous ice in deeper zones.Based on the 1D and 2D velocity profiles, we conclude that generally, the 2-8 m zone contains higher ice connectivity (and therefore less scattered or discontinuous ice), which causes high-velocity zones.This can clearly be seen in 2D plots where a high-velocity zone (and not a continuous layer) exists on two undisturbed permafrost velocity profiles (Figures 4b  and 4c).Therefore, although the temperature can be consistent at higher depths (as shown in Figure 7 using the literature data sets collected near Utqiagvik), the velocity could be different due to changes in ice content and the spatial distribution of ice.While the temperature profiles closely align with each other, there are clear disparities exist in the V s profiles.Such differences might be indicative of pronounced lateral variations in ice content or could reflect changes in the lithology or texture of the underlying sediments.
Ice-wedge polygons occur on nearly all nearshore land surfaces (Kanevskiy et al., 2013) and can be outlined using 2D V s profiles.The formation and degradation of these polygons are linked to climate change, resulting in severe landscape alteration.There are mainly three types of ice-wedge polygons in the tundra between Utqiaġvik and Elson Lagoon: high-centered polygons, flat-centered polygons (incipient polygons), and low-centered polygons.High-centered polygons are shown in Figure 8b, flat-centered polygons in Figure 8c, and low-centered polygons in Figure 8d, surrounded by thermokarst lakes.An early stage of high-centered polygon formation can be seen in Figure 8a.MASW 4 (Figure 8d) is surrounded by coalescent low-center polygons (Lara et al., 2015) and thermokarst lakes based on the satellite view nearby, presenting lower V s compared with high-centered polygon regions (Figure 8b).The landscapes near MASW 4 develop and degrade from flat-centered and high-centered polygons (Nitzbon et al., 2019), showing severe landscape alteration due to climate change.This transformation is referred to as ice-wedge polygon degradation.A water body in the center of the low-centered polygons can change the hydrological regime of polygon nets and lead to the onset of thermokarst activity (Kartoziia, 2019).As shown in Figure 8b, some of the high-centered polygons are developing and connecting, presenting ice-rich permafrost zones with high V s , which have the potential to form thermokarst lakes during permafrost degradation.MASW 2 and 3 cover the polygon centers and troughs of high-centered polygons, while MASW 4 is on the rim between thermokarst lakes.For locations with surface water, V s presents lower values on the top of the permafrost layer than adjacent permafrost.

Influence of Civil Infrastructure on Permafrost
In this section, we discuss the influence of civil infrastructures, including gravel roads and pile foundations, two of the most common civil infrastructures in Northern Alaska, based on seven MASW surveys.Comparison of V s profiles of disturbed permafrost locations (Roadside and NOAA 1) and relatively undisturbed permafrost locations nearby (NOAA 3 and MASW 1) in Figure 4 demonstrates that the ALT is larger in disturbed permafrost due to higher surface temperature.For NOAA 3 and MASW 1, the maximum V s is similar (∼900 m/s), but the highvelocity zone is deeper (∼4 m) in permafrost near civil infrastructure (NOAA 3) compared with MASW 1 (∼1 m).This discrepancy may be due to the disturbed gravel topsoil and higher air temperatures near civil infrastructure, causing diffusive heat transfer from a more absorptive material and resulting in temperature profiles that differ from undisturbed permafrost locations.
The MASW testing location beside a gravel road is shown in Figure 8e (Roadside).In cold regions, dry coarsegrained soil is often used to replace the foundation soil of roadbeds or airport runways to prevent frost heave (Vinson et al., 1996).The gravel fill reduces the road's frost heave and thaw settlement by providing better drainage capability but affects the moisture regime near the gravel road.As shown in Figure 8, high-centered polygons developed near the gravel road, with surface water accumulation next to the road embankment.
Ice-rich permafrost zones can be identified beneath the polygon landscape in the V s profiles shown in Figures 4  and 8.The depth of the ice-rich permafrost zone along the roadside (Figure 4e) is shallower than the nearby tundra location at NOAA 3 (Figure 4g), suggesting the influence of the gravel road.Along the roadside, the ice-rich permafrost zone (∼3 m thickness) is thicker than MASW 1. Different moisture migrations beside the gravel road may cause these differences.In addition to unfrozen water migration as the dominant mode of moisture movement, vapor flux also contributes to frost heaving (Currie, 1983;Farouki, 1981;Smith & Burn, 1987;Teng et al., 2020).Gaseous water (vapor) migrates from the warm and humid side of the soil layer to the cold and dry layer below the closed and impermeable ground surface in coarse-grained soil and then condenses into ice, causing frost heaving (Guthrie et al., 2006;Niu et al., 2017;Zhang et al., 2020).This phenomenon is known as the "pot effect" or "canopy effect".Generally, soil with an initial moisture content of less than 30% is more prone to showing the "pot effect" (Bai et al., 2018).
Pile foundations are the most common building foundation type in Arctic Alaska to overcome differential settlement.Figures 4f and 4g display the 2D V s profiles for MASW surveys under the NOAA building (NOAA 1), and ∼80 m from the building in the tundra (NOAA 3), respectively.As shown in Figure 3, NOAA 1 shows a similar low-high-low V s trend to NOAA 3 (and also MASW 1-4).At depths of 0-2 m, shear wave velocities are slightly different for NOAA 1 and NOAA 3 due to the topsoil of NOAA 1 being gravel, while NOAA 3 is tundra permafrost.At 2-8 m depths, NOAA 3 presents an ice-rich permafrost zone, while NOAA 1 has much smaller V s in this depth range, indicating softer soil.In addition, in NOAA 1, we observed a ∼150 m/s decrease in V s for the ice-rich zone compared to the ice-rich zone at NOAA 3 (farther into the tundra).This lower V s in the ice-rich zone near the building suggests that the pile foundation impacts the soil properties in the surrounding area.Although there is lower V s in the ice-rich zone near the building, the ice-rich zone near the building is more laterally uniform than the ice-rich zone further in the tundra.Because the building had been present in the area for many years, it could have contributed to the thawing and freezing of the surrounding ground, leading to a more uniform distribution of ice-rich soil after years of thermal diffusion of heat from the building.The substantial differences observed at these sites highlight the need to consider the long-term effects of anthropogenic activities on the geological and geotechnical properties of the ground.

Applications in Quantifying Engineering Properties and Designing Infrastructure on Permafrost
Investigating the long-term effect of civil infrastructure on permafrost's stiffness could help improve the engineering design of structures' foundations on permafrost.Permafrost's parameters, such as V S30 and E, are affected by soil types (Coduto, 1999), soil temperature (Ji et al., 2023), and soil's ice content (Fisher et al., 2020).Changes in these soil properties show possible changes in soil types, thermal conditions, and cryostructure.For instance, the stability of the subgrade in permafrost regions, as noted by Anhua (2014), is closely linked to the ice content in the permafrost beneath roadways.Here, we quantitatively analyze V s profiles at NOAA 1 and NOAA 3 locations.
Based on the MASW results, V S30 for locations NOAA 1 and NOAA 3 are equal to 744.2 and 799.5 m/s, respectively.This reduction can affect the soil's ability to carry loads (as the soil at NOAA 1 is less stiff than the NOAA 3 location), leading to greater settlement or deformation under structural loads.Assuming a V p /V s = 1.6 for the site's permafrost layer based on (Ji et al., 2023), and γ and g equal to 19.62 kN/m 3 and 9.81 m/s 2 , respectively, Equations 1-3 result in an elastic modulus of 261.07 MPa at the NOAA 1 location and 301.31MPa at the NOAA 3 location.This indicates a 13.35% reduction in elastic modulus accompanied by an equivalent increase in settlement values, as delineated by Equation 4. While our findings offer useful insights for designing and maintaining infrastructure in polar regions, extending these findings to other permafrost settings should be approached with caution, considering local variations in soil composition, temperature trends, and permafrost degradation.

Conclusions
This study uses 1D and 2D V s profiles from MASW along with temperature measurement, ERT, and permafrost sampling to reveal various permafrost features in Utqiaġvik, Alaska.V s profiles, combined with electrical resistivity models and temperature measurements, can qualitatively characterize active layer, ice-rich permafrost, and cryopeg in permafrost layer, but cannot identify small-scale cryostructures such as ice lenses.V s in the active layer ranged from 240 to 370 m/s (silty peat to silt), while V s in the permafrost layer ranged from 450 to 1,700 m/s (silt to slightly sandy silt) in August 2022.V s profiles demonstrate a consistent vertical low-high-low velocity trend in permafrost.Ice content, ice layers, and ice wedge influence shear wave velocities, with higher V s indicating higher ice content.Low V s permafrost zones may exist across the tundra near Elson Lagoon and east of Utqiaġvik.The V s variation in ice-rich permafrost correlates with ground temperature variation at 0-15 m depths in the study region.This correlation indicates that ice-rich permafrost with higher V s values demonstrates lower temperatures than active layer and ice-poor permafrost.By using ERT, multiple layers can be identified at shallow depths: active layer (200-300 Ωm), cryopeg (10-20 Ωm), and ice-rich permafrost (100-600 Ωm).The presence of ice becomes evident through the analysis of V s and ERT profiles.To strengthen our conclusions, we validate geophysical results with stratigraphy and salinity analyses from permafrost cores.Integrating geophysical, temperature, and core sampling methods offers a reliable approach to evaluating and understanding permafrost spatial variability.
Civil infrastructure can impact permafrost, resulting in a higher active layer thickness and lower V s .The influence of gravel road and pile foundation on permafrost degradation varies.Thicker ice-rich permafrost layers at shallower depths, surface water accumulation, and ice polygon development are identified near the gravel road on permafrost.At the sites with building and pile foundations, lower shear wave velocities are observed at depths shallower than 7 m when compared to nearby undisturbed tundra.The active layer and permafrost are more laterally homogeneous closer to the building compared to nearby undisturbed tundra, and a thinner highvelocity zone exists closer to the building.The resulting V s profiles suggest that weaker ground is present near infrastructure, a factor that should be considered by civil engineers in planning, design, and construction of structures.

Data Availability Statement
• The seismic and ERT data used for geophysical data processing are available through the Arctic Data Center Tourei et al. (2023) https://doi.org/10.18739/A2V40K14Q.• The seismic data were processed using the SeisImagerSW software (GeometricsTM, Version 3.0.),with parameters described in Section 4.1.• The ERT data were processed using Prosys II software (IRIS Instruments) and the Res2Dinv software (Geotomo SoftwareTM, Version 3.59.),with parameters described in Section 4.3.• The temperature data are available through the Arctic Data Center Nicolsky and Wright (2023)

Figure 2 .
Figure 2. The procedure of building V s models from the extracted dispersion curve using the MASW method: (a) The pre-processed shot-gather (Red star represents the shot location), (b) The calculated dispersion image representing Rayleigh wave phase velocity for each frequency (Red dots represent picks at high amplitudes), (c) The extracted dispersion curve from the dispersion image, and (d) V s model inverted from the dispersion data.
02°C.Data records were collected during the field trip in August 2022.Because most of the commonly used construction materials are Journal of Geophysical Research: Earth Surface 10.1029/2023JF007352 TOUREI ET AL.

Figure 5 .
Figure 5. ERT (a) and MASW (b) comparison results at the MASW 3 location.Black lines represent interpreted ERT boundaries based on resistivity values.Zones A-D specify different ranges of electrical resistivity and shear wave velocity values corresponding to various permafrost structures.The white box indicates the borehole sample, and the black box inside the white box shows the transition zone.

Figure 6 .
Figure 6.Stratigraphy of borehole core samples at five locations.

Figure 8 .
Figure 8. Satellite view of MASW testing locations in undisturbed permafrost tundra and disturbed permafrost roadside: (a) MASW 1 in the developing polygon trough area, (b) MASW 2 in the flat-centered polygon area, (c) MASW 3 in the low-centered polygon area, (d) MASW 4 surrounded by thermokarst lakes, and (e) Roadside.