Late Holocene Paleomagnetic Secular Variation in the Chukchi Sea, Arctic Ocean

The geomagnetic field behavior in polar regions remains poorly understood and documented. Although a number of Late Holocene paleomagnetic secular variation (PSV) records exist from marginal settings of the Amerasian Basin in the Arctic Ocean, their age control often relies on a handful of radiocarbon dates to constrain ages over the past 4,200 years. Here we present well‐dated Late Holocene PSV records from two sediment cores recovered from the Chukchi Sea, Arctic Ocean. The records are dated using 26 14C measurements, with local marine reservoir corrections calibrated using tephra layers from the 3.6 cal ka BP Aniakchak eruption in Northern Alaska. These 14C‐based chronologies are extended into the post‐bomb era using caesium‐137 dating, and mercury isochrons. Paleomagnetic measurements and rock magnetic analyses reveal stable characteristic remanent magnetization directions, and a magnetic mineralogy dominated by low‐coercivity minerals. The PSV records conform well to global spherical harmonic field model outputs. Centennial to millennial scale directional features are synchronous between the cores and other Western Arctic records from the area. Due to the robust chronology, these new high‐resolution PSV records provide a valuable contribution to the characterization of geomagnetic field behavior in the Arctic over the past few thousand years, and can aid in developing age models for suitable sediments found in this region.

Over the past decades, a number of studies provided insights into short-term changes of the geomagnetic field in the Canadian and Alaskan Arctic, and as a result, a more complete picture of Holocene paleosecular variation (PSV) in the area has emerged (Barletta et al., 2008(Barletta et al., , 2010Deschamps et al., 2018;Lisé-Pronovost et al., 2009;Lund et al., 2016). Published paleomagnetic secular variation (PSV) records derived from sediment cores from the Chukchi and Beaufort seas appear to display synchronous behavior, although differences up to ∼1 ka exist between some apparently similar features (Lund et al., 2016). The age-depth models of these sediment cores have variable resolution, and are supported by some radiocarbon, tephra, and lead isotope ( 210 Pb) dating, however, for the Late Holocene, they rely on only a few radiocarbon dates. Uncertainties in the geochronologies hinder the development of a PSV master curve for the region that could potentially be used for relative dating of the often microfossil-poor sediment cores recovered from the area.
In order to improve the temporal resolution of paleomagnetic records from this part of the Arctic, Late Holocene PSV records of two sediment cores from the Chukchi Sea are presented. We refine the published age-depth models of the cores by presenting new radiocarbon dates and applying the most recent radiocarbon calibration curve, Marine20 . The chronologies for the past 100-150 years are constrained with downcore mercury and caesium-137 profiles. The new PSV records are compared to geomagnetic field model outputs, and PSV records from sediment cores from the Chukchi and Beaufort Seas.

Regional Setting
The Chukchi Sea is a shallow shelf sea of the Arctic Ocean that occupies the area immediately north of the Bering Strait, a gateway to the Arctic (Figure 1). The opening of the strait at ∼11 cal ka BP  exposed the Chukchi Sea to the direct influence of Pacific waters, which-upon entering the Arctic-split into a number of currents (e.g., Pickart et al., 2016) (Figure 1). As they move northward, these water masses modulate the sedimentation regime of the area. For example, surface sediments from the Chukchi Sea presently are enriched in chlorite and smectite, and have generally high organic carbon (1.5%-2.8% - Stein et al., 2004) and biogenic silica contents (>10% wt), a direct result of inflowing Pacific waters and high seasonal biological productivity (e.g., Swärd et al., 2018). A distinct branch of Pacific waters heads north via the Herald Canyon, an underwater area at the north-western part of the Chukchi Sea ( Figure 1). Sediment cores investigated in this study were collected from the Herald Canyon. The flow dynamics of the canyon exhibit quasi-seasonal variations, with occasional reverse flow (southward) and upwelling of Atlantic waters (Pickart et al., 2010).

Sediment Cores and Their Age Models
The two studied sediment cores (2-PC1 and 4-PC1) were retrieved in 2014 during the SWERUS-C3 (Swedish-Russian-US Arctic Ocean Investigation of Climate-Cryosphere-Carbon Interactions) expedition from relatively shallow waters in the Herald Canyon (Table 1, Figure 1). During the collection of 4-PC1, a 46-cm-long trigger-weight core (4-TWC1) was also obtained. Star-Oddi ® orientation and tilt sensors were attached to the outside of the piston core barrels to monitor tilt through core penetration (O'Regan et al., 2016) (Figure S1a in Supporting Information S1). The cores are composed of gray-olive gray sediments, which have a low and sporadic abundance of calcareous microfossils and mollusks (Seidenstein et al., 2018), and relatively high biogenic silica (∼15% - Swärd et al., 2018) and total organic carbon (TOC) content (0.74%-2.38% - Martens et al., 2019). Detailed descriptions of core sampling, curation, and physical properties measurements were presented by Cronin et al. (2017), Jakobsson et al. (2017), and Pearce et al. (2017). TOC concentration were measured on a Carlo Erba NC2500 elemental analyzer coupled to a Finnigan Delta V Advantage mass spectrometer at the SIL-lab of the Department of Geological Sciences of Stockholm University using freeze-dried and homogenized samples with carbonates removed.
The age-depth model of core 2-PC1, spanning the past 4,200 years, was developed by Pearce et al. (2017), using 17 accelerator mass spectrometry (AMS) radiocarbon datums from mollusks. An ash layer (at ∼7.04 m core depth) from the Aniakchak caldera forming eruption ∼3.6 cal ka BP ago (CFE II) was used to determine local reservoir correction (∆R) values. The coring site was subject to relatively high sedimentation rates of ∼2 m/ka during this period. To further constrain the chronology of the top part of the core, caesium-137 ( 137 Cs) and lead-210 ( 210 Pb) activity levels, and mercury (Hg) concentrations were determined over the upper 0.96 m in 2-PC1. These analyses were performed at the Leibniz Institute for Baltic Sea Research Warnemünde.
Core 4-PC1 is 6.1 m long, and extends to ∼13.5 ka BP or possibly somewhat older. A preliminary chronology was presented by Jakobsson et al. (2017) and Cronin et al. (2017), and was based on 10 AMS radiocarbon datums in total. A hiatus or condensed section is present at 4-4.13 m core depth, manifested as an abrupt lithologic transition. Above the hiatus, three radiocarbon dates constrained the past ∼4,000 years in the core. We obtained an additional 10 radiocarbon ages from fragmented mollusk shells found in the upper 4 m of the core. The caesium-137 profile from the top of the core was used to pinpoint the onset of nuclear bomb testing in the 1950s. Mercury concentrations were obtained every cm from the top 72 cm of the core. The identification of the Aniakchak eruption at 3.73 m core depth (Geels, 2019) was used to establish ∆R for 4-PC1, assuming a constant ∆R over time. In the  Pickart et al. (2016). Inset shows the coring sites of the SWERUS-C3 cores investigated in this study. Basemap: IBCAO version 4.0   Pearce et al. (2017) and Geels (2019) ∆R was determined using the Marine13 calibration curve (Reimer et al., 2013).
Although only separated by 38 km, the depositional settings of 2-PC1 and 4-PC1 are quite different. The sediments of 2-PC1 are mainly hemipelagic in nature, while core 4-PC1 sampled sediments from a drift deposit Swärd et al., 2018). The water depth at the coring site of 4-PC1 (120 m) is approximately twice that of 2-PC1 (57 m), which means that Atlantic sourced waters are sometimes present at the former site and results in different ∆R values for regional marine reservoir age correction Jakobsson et al., 2017;Pearce et al., 2017). The different hydrodynamic conditions present at the two coring sites are also reflected in the sedimentation rates, and it is apparent from the published age-depth models that sedimentation at the coring site of 2-PC1 was almost double that of 4-PC1 over the past ∼4,000 years.
Since publication of the initial age models for these cores, new marine (Marine20) and atmospheric (IntCal20) calibration curves have been introduced Reimer et al., 2020). In the new Marine20 calibration curve the global ocean reservoir correction is about 100-150 years older during the Late Holocene than the previous Marine13 curve (Reimer et al., 2013). Hence, depending on how they were determined, applicable local reservoir corrections will vary considerably when using the Marine20 curve. In this study, previously published radiocarbon ages were recalibrated to produce updated age-depth models. The Marine20 curve is intended for marine samples from non-polar regions only, as variable sea ice and wind strength can complicate air-sea exchange . As a result, the Marine20 calibration curve might not accurately reflect dissolved inorganic carbon ages in polar regions, but this is more likely a problem in older glacial and deglacial time intervals , and potentially less problematic during the Late Holocene. Although the coring sites of 2-PC1 and 4-PC1 are subject to variable summer sea-ice conditions, surface waters in these areas are often free of summer sea ice. To account for the potential limitations in the application of the Marine20 curve, we utilize the IntCal20  curve and the OxCal software of Bronk Ramsey (2009) to calibrate the sample ages, and compare the results with calibrated ages obtained using the Marine20 curve. Prior to each calibration, marine reservoir age offset (∆R) values were estimated by OxCal's built-in Markov chain Monte Carlo (MCMC) analysis using the date of the Aniakchak eruption and 3-4 surrounding radiocarbon dates as was originally done by Pearce et al. (2017) for 2-PC1. As IntCal20 is an atmospheric calibration curve, the resulting ∆R values (hereafter referred to as ∆R_IC) are greater than those from the Marine20 calibration as they include both the reservoir correction (R) and its offset (∆R).
The precise age of the Aniakchak CFE II eruption is not fully established, and Davies et al. (2016) showed that radiocarbon-based ages from Alaska did not agree with those from a Greenland ice core (GRIP) (Abbott & Davies, 2012). Pearce et al. (2017) accounted for the offset between the IntCal13 (Reimer et al., 2013) and GICC05 (Adolphi & Muscheler, 2016) timescales (19 ± 3 years at 3600 BP) to resolve the apparent discrepancy, and proposed 3572 ± 4 cal yr BP for the Aniakchak eruption. The offset between the IntCal20 and the Greenland ice cores chronology is unknown, but it is unlikely to be very different from that of IntCal13, as the underlying Holocene tree ring chronology for IntCal13 and IntCal20 remains the same (Muscheler et al., 2020). Hence, an age of 3572 ± 4 cal yr BP was used in this study for the Aniakchak CFE II eruption.

Paleo-and Rock Magnetic Measurements
To investigate the magnetic fabric of the samples, anisotropy of magnetic susceptibility (AMS) was measured using an AGICO Kappabridge MFK1-FA modular system at the Department of Geological Sciences, Stockholm University. An induced field of 425 Am −1 was used, following calibration with the manufacturer's standard. Data were analyzed with the Anisoft 4.2 data browser (Chadima & Jelinek, 2008).
Paleo-and rock magnetic measurements were conducted on discrete subsamples (plastic cubes with an internal volume of ∼7 cm 3 ) using a 2G Enterprises 760 superconducting rock magnetometer at the Palaeomagnetic Laboratory of Lund University. In total, 320 and 156 discrete subsamples were collected at 2.5 cm increments from core 2-PC1 and 4-PC1 respectively. Measurements included the measurement and demagnetization of the natural remanent magnetization (NRM), anhysteretic remanent (ARM) and isothermal remanent magnetization (IRM), followed by a backfield IRM measurement (all at room temperature). Measurement and subsequent demagnetization of the NRM followed a demagnetization sequence with peak alternating fields (AF) of 5, 10, 15, 20, 25, 30, 40, 50, 60, and 80 mT. The laboratory induced ARM was imparted using a peak AF of 100 mT (with 0.05 mT direct current biasing field), while IRM was imparted with a direct current pulse field of 1 T. The IRM backfield was imparted with a field of 300 mT. Both the ARM and IRM were demagnetized using peak fields of 10, 20, 30, 40, 60 and 80 mT.
Principal component analysis (PCA) of Kirschvink (1980) was used to determine the characteristic remanent magnetization (ChRM) directions (inclination and declination) in both cores. Median destructive field (MDF) of the NRM, and maximum angular deviation (MAD) were also calculated. The MDF is the peak AF required to reduce NRM intensity to half of its initial value. MAD values characterize the precision with which the best fit line, provided by the PCA, is determined (e.g., Butler, 2004). MAD values below 5° could be evidence of well-defined magnetic vectors in Quaternary sediments (Stoner & St-Onge, 2007).
Due to the shallow water depths, and relatively rapid drift of the ship during coring, both cores entered the sediments at an angle that was not perpendicular to the Earth's surface (2-PC1 at 5-18°, and 4-PC1 at 14.5°) (O'Regan et al., 2016), and a correction was applied to counteract this offset. Initially, the core tilt direction was transferred to the sample coordinate system. The tilt azimuth of the core barrel was based on compass readings of the Star-Oddi© sensor and the International Geomagnetic Reference Field (IGRF-12) prediction for the site coordinates. A zero-mean core declination (geocentric axial dipole (GAD) field) was assumed since orientation of the plastic liner relative to the core barrel was unknown. Due to the extensive ship movement during the collection of 2-PC1 and low sampling frequency of the tilt sensor (one sample every 5 s compared to 1 sample per s for 4-PC1), the exact tilt value at the moment of core penetration was difficult to constrain, with values in the range of 5-18° ( Figure S1b in Supporting Information S1). Thus, tilt values of 5° and 18° were assessed to identify which led to realistic tilt correction. The NRM intensity after 30 mT AF demagnetization was normalized by magnetic concentration dependent parameters (ARM 30mT , IRM 30mT , κ), and the ratio was assessed to identify an appropriate normalizer capable of providing an estimate of relative paleointensity (RPI). Different AF demagnetization steps (20-60 mT range) yielded the same result.
The results of radiocarbon, and caesium-137 dating are summarized in Tables 3a and 3b and Figure 2. The activity of 210 Pb unsupported was very low, and hence it was not used as a chronometric marker. Caesium-137 entered sediments from global fallout from atmospheric nuclear weapons testing in the 1950s and early 1960s, with major global inputs in 1959 and 1963/64 (e.g., Ritchie & McHenry, 1990). Hence, a calibrated date of −9 ± 5 cal yr BP (equivalent to 1959 CE) was used in the age modeling of the cores. Similar increases in mercury concentration observed near the core tops in both records ( Figure 2) were assumed to be synchronous, and were set as stratigraphic tie points between the records. The Aniakchak CFE II tephra provided a further tie point. For the OxCal Bayesian age-depth model, a P-sequence routine (Bronk Ramsey, 2008;Bronk Ramsey & Lee, 2013), which assumes deposition to be Poisson distributed, was employed ( Figure 3). Two outlier radiocarbon ages were removed (Beta-487861 and Beta-512060) from the age model of core 4-PC1, and a previously published radiocarbon age (LuS11278, Jakobsson et al., 2017) was also omitted as it would return a "post-bomb" carbon date once the reservoir effect was removed.
Recalibration of radiocarbon dates in core 2-PC1 with the Marine20 and IntCal20 calibration curves resulted in mean ages statistically indistinguishable from the previously published ages  for all samples Note. Also shown are previously published estimates by Pearce et al. (2017) and Geels (2019) based on the Marine13 curve.
a Includes global reservoir age (R) and its regional correction (∆R).

Note.
Radiocarbon ages were calibrated with the Marine20  and IntCal20  calibration curves. Also shown are previously published mean Marine13 (Reimer et al., 2013) calibrated ages and uncertainties from Pearce et al. (2017). Aniakchak CFE II tephra horizon is present at 7.040 ± 0.075 m . pMC = percent modern carbon. Samples in gray italics were not used in the age modeling.  (Tables 3a and 3b). As the calibrations resulted in similar age-depth models, hereafter the age model based on the Marine20 calibration (Figure 3) was used.

Magnetic Fabrics
Principal axes (K max , K int , K min ) of the susceptibility ellipsoids are shown (without tilt correction) in Figure 4. K max and K int are generally dispersed in the horizontal plane in both cores, but K max declinations exhibit limited preferential diagonal direction (∼135-315°) in the top sections (Sections 1), and in the trigger-weight core (4-TWC1). K min axes of subsamples from all but the top section of core 2-PC1 cluster and trend toward the vertical, characteristic of horizontally bedded sediments, but in Section 1 of 2-PC1, and in cores 4-PC1 and 4-TWC1, they shallow and are more dispersed. Mean K min inclination of 2-PC1 is consistent with a core tilt of 5°, however as previously discussed, correcting for a core tilt of this magnitude results in relatively shallow palaeomagnetic directions. The relationship between the degree of anisotropy (Pj) and the shape parameter (T) indicates that the majority of subsamples have an oblate ellipsoid (more so in 2-PC1), typically attributed to disturbance-free horizontal bedding planes.

Paleo-and Sediment Magnetism Results
Two separate components of the NRM could be identified in the cores (including 4-TWC1): a low coercivity component was readily removed by AF of 10-20 mT, which revealed a second stable component with constant NRM directions between 20 and 80 mT ( Figure S3 in Supporting Information S1). ChRM directions were based on the results from AF demagnetization steps between 20 and 80 mT. In both cores, MAD values are generally below 5° (74% and 92% of directions), which is indicative of well-defined characteristic components (e.g., Stoner & St-Onge, 2007), assuming the lack of sampling induced deformation of the magnetic fabric.
The GAD theorem predicts a field inclination of 81° for the coring sites. Due to the uncertainty in the extent of core tilting of 2-PC1, two different core tilt values (5° and 18°) were applied during tilt correction. These resulted in mean ChRM inclination values of 63° ± 11° and 73° ± 10° respectively ( Figure S4 in Supporting Information S1) for 2-PC1. Although the AMS data showed signs of a tilt consistent with 5°, from a paleomagnetic perspective a correction leading to higher mean inclination value is more likely, hence we applied a correction for a core tilt of 18° for core 2-PC1. For 4-PC1, a core tilt of 14.5° was used ( Figure S4 in Supporting Information S1) and this resulted in mean ChRM inclination of 65° ± 17°. The tilt-corrected ChRM inclination records of both cores are somewhat shallower than the GAD prediction ( Figure 5). Section 1 of core 4-PC1 is characterized by particularly shallow values (37.6° ± 7.5°), and these low values can also be observed in the ChRM inclination profile of the corresponding trigger-weight core, 4-TWC1 ( Figure S5 in Supporting Information S1). The fact that both the piston core and trigger-weight core show the same shallower than predicted inclinations suggests that this is unlikely to be a result of vertical stretching of sediment or other coring induced deformations, and is more related to primary depositional processes (current induced) or post-coring deformation of relatively soft sediments. This is further confirmed by ChRM declination values (following setting core mean to zero) from the uppermost core sections, which are distinctly shifted from the rest ( Figure 5).
To improve the statistical quality of the paleomagnetic directions, the recommendation of Deenen et al. (2011) was followed and ChRM values were converted to a virtual geomagnetic pole using the paleomagnetism. org platform of Koymans et al. (2016). A 45° pole cut-off was then applied to reject data, which resulted in the rejection of 20 subsamples (equivalent to 6.5%) in core 2-PC1, and 24 subsamples (15.8%) in 4-PC1.
Mean values of the MDF NRM (30 ± 14 mT and 40 ± 8 mT) in the two cores indicate the dominance of low-coercivity minerals. Both records are characterized by a distinct interval of reduced MDF NRM values exhibiting high variability ( Figure 5). The onset of this interval is synchronous in the cores and corresponds to ∼3.1 ka, but the interval lasts longer in core 2-PC1 (∼1.7 ka cal yr BP) than in 4-PC1 (∼2.4 ka cal yr BP). The presence of these intervals could signal diagenetic alteration or dissolution of the magnetic minerals.
Concentration dependent parameters, such as low-field susceptibility (κ), ARM and IRM largely covary, and their profiles are also characterized by a distinct interval with reduced values in both cores ( Figure 6). This interval is also present in the profiles of the grain size sensitive parameter, ARM/κ, indicating increased magnetic grain size at this depth. The range of variability, however, remains within one order of magnitude for all parameters, and is largest for the ARM in both cores. The S-ratio (in this study defined as the ratio of the backfield IRM 300mT and the IRM 1T ) exhibits minimal changes in both cores, and its mean values of −0.96 ± 0.01 (2-PC1) and −0.95 ± 0.01 (4-PC1) are characteristic of low-coercivity ferrimagnetic minerals (e.g., Thompson & Oldfield, 1986).

Relative Paleointensity
RPI estimates can approximate magnetic field intensity variations, as long as the RPI data are characterized by high magnetic stability of the natural remanence, moderate (less than an order of magnitude) changes in   with ∆R = 330 ± 41 years (2-PC1) and ∆R = 206 ± 67 years (4-PC1). Constant sedimentation rate was assumed for the intervals between the core bases and overlying radiocarbon dates. Laboratory identifiers of previously published radiocarbon ages (obtained from Cronin et al., 2017;Jakobsson et al., 2017;Pearce et al., 2017) are shown in gray italics, and outliers are displayed in red. Dark and light shading (blue/green) represents 1σ and 2σ age ranges. concentration, magnetic mineralogy and grain size, and are supported by a robust independent chronology (e.g., Tauxe & Yamazaki, 2015). Normalizing the NRM 30mT by IRM 30mT and κ produced RPI estimates that resembled bulk rock-magnetic parameters (ARM, IRM and κ), but, normalization of the NRM 30mT by ARM 30mT produced a profile that is somewhat different from these (Figures 6 and 7). However, it is notable that intervals of high NRM 30mT /ARM 30mT (around ∼2000-3000 cal yr BP) are associated with low magnetic concentration and Figure 5. Characteristic remanent magnetization (ChRM) declination and inclination, natural remanent magnetization (NRM) intensity after 30 mT demagnetization, median destructive field of the NRM, maximum angular deviation of the ChRM fits, sediment mean grain size, and wet bulk density (ρ) and total organic carbon content (%) for cores 2-PC1 and 4-PC1 (prior to data rejection). Grain size data from Swärd et al. (2018). Red dashed lines on the ChRM inclination plots represent 81°, the expected inclination at the coring site for a geocentric axial dipole. Dashed horizontal lines mark core sections.
coarser magnetic grain size, signaling a clear environmental influence. Therefore, a valid RPI record could not be produced for the cores.

Age Model Refinement via Automated Stratigraphic Alignment
Considering the geographical proximity of the two coring sites (∼38 km apart), highly similar geomagnetic field variation, thus PSV records, can be anticipated. Diagenesis potentially can alter these, and the presence of an interval with oscillating MDF NRM values could suggest that it might also be acting on the PSV records of the two cores. Variable MDF NRM values were previously observed in sediments from the Chukchi-Alaskan margin by  Brachfeld et al. (2009), who linked their oscillation to changes in sea ice cover. They proposed that intervals with low MDF NRM values are characterized by reduced amounts of greigite (and reduced sea ice cover), and hence, an NRM dominated by magnetite capable of accurately recording palaeomagnetic directions. Intervals with low(er) and oscillating MDF NRM values in the cores could reflect the increased influence of magnetite, but could also result from variable concentration and grain size of magnetic minerals. Although, the differences in sedimentation rates between the cores are reflected in the lower resolution magnetic record of core 4-PC1, nevertheless, the ChRM inclination records of the two cores display good agreement, which allows their synchronization (Figures 5 and 8a). While the ChRM declination profiles of the cores also exhibit similar variability, declination values near the geomagnetic pole are generally more prone to significant fluctuations, thus here we focus on the inclination records only. The stratigraphic alignment of 2-PC1 and 4-PC1 was modeled using an automated algorithm for probabilistic inversion. The inverse problem is formulated using a Bayesian framework that samples the full range of possible alignment scenarios between a given set of input and target proxy records and explicitly builds in prior chronostratigraphic information. The numerical approach builds on previous work using a hidden MCMC model for automated synchronization of proxy data, which was presented in a seminal paper (Muschitiello et al., 2020) and has been successfully applied on a variety of paleoceanographic records Sessford et al., 2019;West et al., 2019).
In this study we performed a multi-parameter alignment that simultaneously correlates the stratigraphy of 2-PC1 and 4-PC1 using two independent proxy data sets. We used the ChRM inclination records and the AMS 14 C series from core 4-PC1 as inputs and their counterpart data from core 2-PC1 as targets (Figures 8a and 8b). In addition, the model considers the stratigraphic tie points provided by the Aniakchak tephra and by the 137 Cs and Hg profiles, as well as their uncertainty in the depth domain (Tables 3a and 3b).
Finally, the method weighs the probability of any given alignment based on the misfit between the inputs and the targets, and ultimately "finds" a sample of alignments of 4-PC1 onto the 2-PC1 stratigraphy that are consistent with the common chronostratigraphic markers. The algorithm hinges on the assumption of direct synchrony of fluctuations in the ChRM inclination signals and AMS 14 C ages at both coring sites and circumvents the limitations generally associated with subjective and point-wise visual alignments, thus providing a reproducible and continuous alignment that accounts for potential uneven compaction/expansion in sediment cores.
The simulation of the final stratigraphic alignment was run for 10 6 iterations after discarding a burn-in time of 10 5 MCMC steps. The sample from the remaining iterations was used to estimate the posterior median alignment The ChRM data are presented in normalized units. Prior to synchronization, the ChRM inclination data were high-pass filtered using a cut-off period of 10 cm to remove high-frequency noise, and the input 14 C dates of 4-PC1 were corrected by applying a constant ∆R offset of 125 years between the two coring sites. Unfiltered data are also shown. (c) Posterior median (red line), weighed mean (blue line), and pointwise 95% credible intervals (shading and blue dashed lines) of the depth-depth relationship between 4-PC1 and 2-PC1. Stratigraphic markers used to constrain the automated synchronization are also shown. Note that each stratigraphic tie point is associated with an error in the depth domain that is incorporated into the alignment model. that is, the Aniakchak tephra horizon was assigned an uncertainty of ±7.5 cm (2σ) for 2-PC1  and ±5 cm for 4-PC1, whereas the 137 Cs and Hg markers were assigned a conservative uncertainty of ±5 cm for both cores. The largest alignment error (∼10-20 cm) is observed in the interval ∼3-5 m (on the target 2-PC1 depth scale) that is where 4-PC1 lacks 14 C constraints, thus making the synchronization relatively more uncertain.
between 2-PC1 and 4-PC1 and its credibility bands, which reflect the uncertainty of the overall alignment ( Figure 8c).
The overall agreement between the inclination ( Figure 8a) and declination (not shown) profiles and radiocarbon dates (Figure 8b) implies that it is most likely that "true" geomagnetic field behavior was captured in both records despite sedimentological differences, and the potential impact of diagenetic alteration.

Comparison With Geomagnetic Field Model Outputs
Following synchronization of the ChRM inclination records of the cores, the inclination, and declination profiles were compared to geomagnetic field model predictions. Two models were selected: COV-ARCH (Hellio & Gillet, 2018), which utilizes volcanic and archaeomagnetic data only, and "pfm9k.1a" (Nilsson et al., 2014), which additionally incorporates sedimentary data in the model inputs. Prior to comparison, the inclination and declination records of both cores were smoothed using locally weighted scatterplot smoothing (with a 5-point span). ChRM inclinations from both cores generally conform very well to field model outputs, and a reasonably good match can be observed between the inclination profiles and the models for the 800-2600 cal year BP period ( Figure 9a). A number of prominent inclination minima can be observed, particularly at around 1150, 2580 and 3060 cal yr BP.
Despite an overall good match of the ChRM inclination profiles to geomagnetic field models, the mean ChRM inclination of 4-PC1 is ∼10° shallower than that of 2-PC1, and is almost 20° lower than the 81° predicted by the GAD theorem for the coring sites. As the cores were in close-proximity, the influence of disparate non-dipole field components is unlikely to explain such a difference between the two cores. Swärd et al. (2018) showed that current-related sorting and variable bottom-current strength were evident in core 4-PC1 -consistent with its recovery from a drift deposit on the flank of Herald Canyon , hence the shallow inclinations may in part be attributed to bottom water current activity, unless the NRM is primarily a post depositional remanent magnetization. The influence of currents is seen in the results of the anisotropy of magnetic susceptibility analyses (Figure 4), which show that K min axes of subsamples from cores 4-PC1 and 4-TWC1 are dispersed and generally off vertical, and inclinations of K max and K int axes are mostly further away from 0° than in core 2-PC1 (even when corrected for core tilt). These observations indicate that directions of the principal AMS axes deviate from those expected for horizontally bedded sediments. Although coring and subsampling processes can also introduce such biases (Snowball et al., 2019), these are unlikely, as the magnetic fabric of core 2-PC1, apart from Section 1, lacks evidence for such disturbance, and we have a priori information that sediments in 4-PC were affected by variable bottom water currents (Swärd et al., 2018). The presence of disturbance in Section 1 of 2-PC1, and the particularly shallow ChRM inclination values from Section 1 (51 cm) of 4-PC1 and 4-TWC1 Figure 9. Comparison of the Characteristic remanent magnetization (ChRM) (a) inclination, and (b) declination profiles of cores 2-PC1 and 4-PC1 to geomagnetic field model outputs for the core locations as predicted by the pfm9k.1a (Nilsson et al., 2014) and COV-ARCH (Hellio & Gillet, 2018) field models. The records are projected on an age scale (following synchronization). Locally weighted scatterplot smoothing (with a 5-point span) was applied to the ChRM inclination and declination profiles. Uncertainty estimates (red and gray shading) for the inclination and declination profiles were calculated by converting maximum angular deviation values to ⍺95 (Khokhlov & Hulot, 2016), which were subsequently converted to standard deviation of inclination and declination (Suttie & Nilsson, 2019). Blue shading marks 1σ uncertainty associated with COV-ARCH.
could reflect the higher water content and less consolidated nature of sediments at the tops of these cores, which may have been affected by storage/sampling rather than being artifacts of the coring process itself. These intervals cover the past ∼250 years, implying that there remains considerable uncertainty in the PSV records over this period.
ChRM declination values from near the core tops also clearly diverge from model outputs in both cores (Figure 9b), which might be tightly linked to the disturbances observed in the magnetic fabric as discussed above. Paleomagnetic declination data from the polar regions often display large variability as inclination values approach 90°, with uncertainties present in both the recording and sampling/sub-sampling process, and this can sometimes make their comparison to model outputs less successful. Nevertheless, the deeper lying ChRM declination profile of the cores appear to have captured the easterly swing in declination around 2200 cal yrs BP that is predicted by both models.

Comparison With Other Western and Alaskan Arctic Records
Published paleosecular variation records from the Chukchi and Beaufort Seas (e.g., Barletta et al., 2008Barletta et al., , 2010Darby et al., 2012;Deschamps et al., 2018;Lisé-Pronovost et al., 2009;Lund et al., 2016) indicate some synchronous changes of the geomagnetic field over the past ∼10 kyr. Common features identified by earlier studies (e.g., Barletta et al., 2008;or Figure 9 in Deschamps et al., 2018), can also be observed in the ChRM inclination profiles of 2-PC1 (and 4-PC1) ( Figure 10). For example, the inclination peak at ∼1850 cal yr BP or the inclination minimum at ∼2600 cal yr BP can be correlated to those in cores from both the Eastern Chukchi Sea (e.g., core HLY0501-6JPC) and the Beaufort Sea (core 2004-804-803), suggesting similar field behavior over relatively large regional scales (>1000 km). There are slight discrepancies between published ages for these features (e.g., Table 3 in Deschamps et al., 2018), however these are not entirely unexpected. Phase offsets between similar PSV features could be linked to the general westward drift of the geomagnetic field at high latitudes (Nilsson et al., 2020), potential differences in paleomagnetic lock-in-depths, the scarcity of absolute ages available in the other cores for the 0-4,500 years BP period, and the variability in the utilized local marine reservoir corrections (which range from zero to several hundred years). Not only is it difficult to determine the exact reservoir correction to be used in Arctic sediments, but ∆R values can also change over time, and recent studies on sediment cores from the area highlighted the complexity of the issue (e.g., Darby et al., 2012;Keigwin et al., 2018).
The apparent synchronicity of PSV features can, however, help resolve some of the uncertainties in defining appropriate local reservoir corrections, and help identify changes in sedimentation rates that may not be captured in many of the other less well dated records. For example, inclination features from cores 2-PC1 can be shown to align with those from core 2004-804-803 (hereafter core 803) collected from the Beaufort Sea (Figure 1). Prior to correlation, ChRM values were converted to virtual geomagnetic poles (Deenen et al., 2011), and a 45° cut-off was applied to reject data that were likely outside the range of reasonable secular variation. (b) Age-depth models for core 2004-804-803: original age-depth model based on Barletta et al. (2008) with ages calibrated with the Marine04 curve (Hughen et al., 2004), and age-depth model based on correlation with core 2-PC1 using the Marine20 curve . Asterisks and circles mark radiocarbon dates, and inclination tie points respectively.
Following tie point correlation of the two records (Figure 11a), ∆R values can be estimated for core 803. By recalibrating the 4 radiocarbon ages with the Marine20 curve constraining the Late Holocene in core 803 (Barletta et al., 2008), and applying OxCal's Bayesian modeling (Bronk Ramsey, 2009), an estimated ∆R value of 283 ± 60 years is obtained for core 803. This estimate and the resulting age-depth model compare well with that of Barletta et al. (2008) (Figure 11b), who applied a regional reservoir correction of 400 years for calibration of the 14 C dates using the Marine04 curve (Hughen et al., 2004).

Conclusions
The new paleosecular variation records from the SWERUS-C3 cores, 2-PC1 and 4-PC1, provide well-dated palaeomagnetic records from the Herald Canyon in the Chukchi Sea for the Late Holocene. Palaeomagnetic directional data from the cores show similar variability, allowing magnetic age refinement of the core with relatively lower sedimentation rate (4-PC1). The data are largely consistent with time-varying geomagnetic field models, such as pfm9k.1a (Nilsson et al., 2014) and COV-ARCH (Hellio & Gillet, 2018), but also show evidence for larger amplitude variations in inclination that the current range of models are unable to capture due to limitations in chronologies and limited number of paleomagnetic records in the polar regions. Our results provide important information on geomagnetic field at high northern latitudes and may shed light on recent observation of persistent westward drift at high latitudes of the core (Nilsson et al., 2020). Grounded by robust age-depth models these PSV records could underpin a regional Holocene reference curve for investigating geomagnetic field behavior at high latitudes, and for the palaeomagnetic age refinement of sediment cores lacking reliable age constraints for the Late Holocene. The study also highlights the importance of recovering and studying multiple cores for paleomagnetic studies, especially at high latitudes, to overcome the well-known issues associated with coring/handling of the cores, disturbance in high water content, diagenesis, and currents.