A Millennial‐Scale Oscillation in Latitudinal Temperature Gradients Along the Western North Atlantic During the Mid‐Holocene

Changes in vegetation in North America indicate Holocene shifts in the latitudinal temperature gradient along the western margin of the North Atlantic. Tree taxa such as oak (Quercus) and birch (Betula) experienced opposing directions of change across different latitudes consistent with changes in temperature gradient steepness. Pollen‐inferred temperatures from 34 sites quantify the gradient changes and reconstruct a long‐term northward steepening in summer and southward steepening in winter. From 4.8 to 3.8 ka, an oscillation in tree distributions interrupted the long‐term trends as a steep temperature gradient developed north of 43.5°N. The shift likely limited cold outbreaks to the south, producing anomalously high summer temperatures at 42–43.5°N, and enabling a northward expansion of oak forests. The forest and temperature gradient changes appear consistent with orbital and ice sheet forcing as well as millennial variability in the North Atlantic pressure field analogous to the North Atlantic Oscillation on interannual time scales.

potentially because shifts in gradients and frontal movements can create opposing changes and areas of no change in close geographic proximity.
The sentivity of Holocene tree distributions to both growing season warmth and winter freezing temperatures can help reconstruct the past temperature gradients (Prentice et al., 1991;Woodward, 1987).In eastern North America near the North Atlantic, an ecotone between mixed and broadleaf forests expressed a close association with temperature at the time of European colonization (Cogbill et al., 2002;Thompson et al., 2013).It angled across the northeast U.S., approximating the orientation of today's 8°C mean annual isotherm with respect to latitude, topography, and proximity to the ocean (Cogbill et al., 2002).However, its position shifted throughout the Holocene (Oswald et al., 2018) and the earliest palynological studies inferred temperature changes from these ecotone dynamics (Deevey, 1939).Because direct comparisons reveal a tight coupling between fossil pollen and independent paleoclimate data (Shuman et al., 2019), the dense network of fossil pollen stratigraphies in the northeast U.S. (Figure 1a) may help to document patterns and timescales of Holocene changes in the latitudinal temperature gradient (LTG) where it is sensitive to the NAO (Figures 1b-1e).
Oak (Quercus spp.) and hickory (Carya spp.) tree populations dominate the broadleaf forests south of the ecotone and may be particularly useful indicators of past changes in the LTG.The tree distributions define the ecotone where they trade-off in dominance with cool-tolerant taxa such as birch (Betula spp.) and hemlock (Tsuga canadensis) to the north (Cogbill et al., 2002).Additionally, oak and hickory differ in their seasonal temperature sensitivities.They have similar distributions in the northeast U.S. today, but strikingly different historic distributions in the central U.S. where continental summer and winter temperature gradients differ from each other; cold continental winters limit hickory's northern extent more than that of oak (Paciorek et al., 2016(Paciorek et al., , 2021;;Prentice et al., 1991).Changes in the geography of these trees therefore may provide insight into Holocene changes in the LTG in both summer and winter, including in association with a long-debated regional decline in hemlock abundance (Booth et al., 2012;M. B. Davis, 1981;Foster et al., 2006).
Here, we use a dense network of 34 pollen stratigraphies to infer summer and winter temperature changes at different latitudes across the northeast U.S. Mean annual temperatures range from <2°C to >12°C across the study area today (Figure 1a), but we reconstruct changes both in the absolute range and the shape of the gradient.We evaluate the reconstructions by comparing them with the regional histories of oak, hickory, birch, and hemlock.The major patterns of change highlight the ecological significance of millennial-scale climate variability near the North Atlantic.

Methods
Detailed fossil pollen records from Connecticut, Massachusetts, New Hampshire, New York, Pennsylvania, and Vermont were obtained from the Neotoma Paleoecological Database (Williams et al., 2018).The records were selected from areas that span the full north-south temperature contrast between 40 and 45°N latitude.This range of latitudes was possible from 70 to 77°W longitude, but not in places to the east, such as Maine.To ensure a wide sampling of the regional patterns (particularly in the north), we included all radiocarbon dated records that span at least the past 8 ka (Table S1 in Supporting Information S1).We account for uncertainities tied to varying chronological control and sample frequency, and reduce local ecological noise, through as much replication as possible within latitudinal zones and averaging the different time series.
Temporal resolution of the records ranges from >750 to <10 years per sample (95% range) with a median resolution of 124 years per sample.Chronologies were updated using Intcal20 in bchron (Parnell et al., 2008;Reimer et al., 2020).Past comparisons of radiocarbon ages from bulk sediment and terrestrial macrofossils in the same cores in this region reveal few differences, reducing the need to remove older records from the synthesis (Marsicek et al., 2013).Although no records stand out as substantial outliers, we identify one highly-resolved record from each of three latitude zones to verify that the major LTG patterns do not depend upon confounding elevation-temperature relationships, sampling resolution, or age control.The three highly-resolved records have similar elevations (300-400 m), sampling frequency (1σ range of 45-158 years/sample) and age control measured by accelerator mass spectrometry; they derive from Sutherland Pond, New York (Maenza-Gmelch, 1997a), Little Pond, Royalston, Massachusetts (Oswald et al., 2018), and Knob Hill Pond, Vermont (Oswald et al., 2018).
We reconstructed mean summer (June-August) and winter (December-February) temperatures using the modern analog technique following Marsicek et al. (2013).We compared each fossil sample to all possible modern Writing -review & editing: Ioana C. Stefanescu, Laurie D. Grigg, David R. Foster, W. Wyatt Oswald analogs from North America east of 95°W using the squared-chord distance metric (Overpeck et al., 1985).The modern pollen samples and their associated temperatures derive from Whitmore et al. (2005).For each fossil pollen sample, the temperatures of the best seven modern analogs were averaged to produce the reconstructed temperature using the R package, rioja (Juggins, 2019).
Hierarchical cluster analysis of the summer temperature reconstructions (hclust based on Euclidean distances in R) (R Core Development Team, 2022) was used to group the individual records into three ensembles.The clusters were determined from the consistency in their summer temperature reconstructions, which break out by latitude with northern (>43.5°N),central (42-43.5°N),and southern (<42°N) groups of sites, except that eastern coastal sites cluster with the southern group (Figure 1a).The grouping is consistent with the angled orientation of the LTG today because of its interaction with both topography and proximity to the coast (Figure 1a).The LTG reconstructions were based on the averages of all individual records within each group, but we confirm the patterns using the three records with similar elevations of 300-400 m.Reconstructions were first interpolated to 50-year intervals to enable averaging.
The position of the steepest portion of the LTG was calculated as the change in temperature between latitude bands by measuring the mean temperature differences among clusters.We calculated both the mean central-northern (C-N) and southern-central (S-C) temperature difference through time.The two segments were then compared to detect where the LTG was steepest by subtracting the C-N and S-C differences from each other.Positive values in the C-N minus S-C difference indicate a steeper gradient north of 42.75°N, the mid-point of the region, and negative values indicate a steeper gradient to the south.

Reconstructed Temperature Trends
All of the temperature reconstructions indicate significant early Holocene warming in both summer (Figure 2) and winter (Figure 3).At ca. 11.7 ka, the range of summer temperatures increased from 15-16°C to 17.5-18.5°Cwith all latitudes warming by similar amounts.Summer warming peaked in the central and northern sites by 8 ka  75°N (1960, 1977, 1988, 1999, 2001) and (c) north of 42. 75°N (1954, 1959, 1968, 1991, 2007).Lower maps represent (d) the associated composite mean temperature differences (LTG steep in north minus steep in south) and (e) the correlation coefficeient, r, between the North Atlantic Oscillation and mean annual temperature (T) in each U.S. climate division (Vose et al., 2014).Modern mean annual temperatures across the region (a) are based on 1991-2020 PRISM normals (PRISM, 2020); symbol colors indicate the geographical cluster containing each fossil site.See Supporting Information S1 for composite anomaly methodology.at 18-19°C, but southern sites continued to warm, reaching 21°C at 5.5 ka (Figure 2).In winter, all three regions continued to warm until 5.5 ka (Figure 3).After 10 ka, however, the southern sites warmed faster in winter than the central and northern sites; by 5.5 ka, southern sites reached near freezing (0°C) even though northern areas only warmed to −5°C (Figure 3).
After the mid-Holocene, summer temperatures declined in both the northern and southern regions, but not at central sites, which were affected by a significant millennial-scale variation at 4.8-3.8ka (Figure 2).Initially, at 6.2 ka, northern sites cooled by 0.5°C in both summer (Figure 2) and winter (Figure 3).Then, by 4.8 ka, central sites warmed 0.5°C in summer (Figure 2) and 1°C in winter (Figure 3).Some northern sites then warmed again by a similar amount from 3.4 to 1.8 ka (thin blue lines, Figures 2c and 3c), but rapid summer cooling of 0.5°C at 1.8 ka further prevented most northern sites from returning to their previous high temperatures.The three best-resolved records from similar elevations (thin lines, Figures 2d and 3d) show patterns consistent with the inherently smoothed regional means (bold lines, Figures 2d and 3d).

Reconstructed Gradient Changes
Qualitative interpretation of the fossil pollen record supports the infered gradient changes.Oak, which flourishes today where summer temperatures exceed 17°C (Williams et al., 2006), reached its maximum in the north by 9-8 ka when birch, a dominant northern pollen type, declined to a minimum (Figure 2, green lines).Northern oak abundance then fell to near zero, even as it continued to increase in the south where birch was uncommon.The opposite direction of oak and birch change in the north and south from 9 to 5.5 ka, observed over many sites, indicate a widening difference in growing conditions.The temperature reconstructions based on the full pollen assemblages suggest that the overall temperature difference between north and south increased by 1-1.5°C (Figure 2e).with the mean differences between clusters (e) and the north minus south difference in gradient steepness ("Gradient Diff." as the black line in (f).An inverted June insolation curve for 30°N (Berger, 1978) is normalized to the gradient difference scale for comparison (thin line, f).Vertical gray lines mark major time periods of change with the bold line at 5.5 ka indicating the median beginning of the regional hemlock decline (see Figure 3).The vertical gray line at 4.8 ka is unlabeled.
The history of the central region differs from both north and south.Oak abundance also reached an early Holocene peak there as the ecotone shifted far north, but oak then declined before forming a second millennial-scale peak during a northward shift in the ecotone from 4.8 to 3.8 ka (Figure 2b, green lines).The distinctive peak developed when oak abundance remained low in the north.The vegetation gradient thus steepened to the north because birch, rather than oak, reached a mid-Holocene maximum (Figure 2c).
Unlike the oak and birch changes, the regional hemlock decline affected the central and northern areas similarly (Figures 3b and 3c).Average birch abundance initially increased in both areas after 5.5 ka, but it then declined in the central region as oak increased by 4.8 ka (Figure 2b).Beech (Fagus grandifolia), another cool-tolerant tree, experienced a similar brief peak from 5.5 to 4.8 ka at central sites (Figure 4a).The mid-Holocene warmth that accompanied the subsequent oak maximum reduced the summer temperature gradient in the south (S-C) by 0.5°C but increased the northern gradient (C-N) by >1°C (Figure 2e).
The reconstructed gradients also changed in winter.The patterns appear qualitatively in the Holocene history of hickory (Carya), a genus of broadleaved deciduous trees usually restricted to areas today with mean winter temperatures greater than −4°C (Williams et al., 2006).Unlike oak, hickory only became important in the study area after ca. 8 ka (green lines, Figure 3).Its history parallels the slower long-term increase in temperatures during winter compared to summer.As oak pollen and summer temperatures declined in the late Holocene, hickory pollen retained its earlier abundance suggesting stable winter temperatures.
The increase in hickory in the south after 8 ka, but not in the north, marks a widening latitudinal vegetative difference.The winter temperature reconstructions generated from the full pollen assemblages show a persistent C-N difference of ∼2°C (thin blue line, Figure 3d), but a widening S-C difference from near-zero to >4°C (thin red line, Figure 3d).An exception to this pattern developed from 4.8 to 3.8 ka when hickory abundance increased in the central region as oak-hickory assemblages expanded into central highland sites.As with oak, the hickory increase followed the regional decline in hemlock by 700 years on average (Figure 3b) and did not develop in the north where the reconstructed winter temperatures declined (Figure 3c).The winter reconstructions indicate

Changing Latitude of Steepest Gradients
The difference in the steepness of the two summer temperature gradients (C-N minus S-C, Figure 2f) indicates: • A southward shift in the steepest part of the gradient before 6 ka as southern areas warmed, but northern areas cooled; • A northward shift from 4.8 to 3.8 ka when the northern (C-N) gradient was steeper than the southern (S-C) gradient by up to 0.8°C, similar to the cool, late-Holocene pattern since 1.8 ka (see Figure S1 in Supporting Information S1 for a geographic representation of the temperature gradient); and • A nearly uniform regional gradient (a near zero difference, Figure 2f) between these two periods from 3.8 to 1.8 ka.
In winter, because the northern (C-N) gradient remained near 2°C and the southern (S-C) gradient increased from 0°C to >4°C, the steepest gradient was initially in the north, but shifted to the south by the late-Holocene.The northward oscillation in the position of the steepest segment interupted the trend from 4.8 to 3.8 ka when C-N was steeper than S-C by up to 1.6°C (Figure 3f).

Interpreting the Long-Term Changes
At the time of European colonization of the northeast U.S., the LTG determined a sharp contrast between northern mixed forests and oak-hickory forests to the south (Cogbill et al., 2002).The reconstructed shifts in the LTG, therefore, provide a plausible climatic explanation for latitudinal differences in the regional vegetation history.A weak north-south gradient in summer before ca.8.2 ka (Figure 2e) helps to explain the early peak in oak abundance at northern sites (green lines, Figure 2c); the polar front was likely furthest north at the time and similar summer temperatures across the region enabled oak to spread widely.When the summer LTG then steepened (Figure 2e), oak readily declined in the north where sites began cooling (Figure 2c) while also expanding in the south (Figure 2a).Tree population increases south of a non-advancing range limit are not uncommon (Lloyd, 2005), and should be expected if the LTG steepened in place of uniform temperature increases or declines regionwide.Birch followed the opposite pattern as oak until the mid-Holocene, widening the vegetation difference (Figures 2a-2c) as the LTG steepened across all latitiudes.
Different directions of insolation forcing explain a northern steepening of the LTG in summer (Figure 2f), but a southern steepening in winter (Figure 3f).Early-to-mid Holocene insolation anomalies acted to weaken the summer LTG by preferentially heating northern continental areas, especially compared to coastal areas in the southeast.The reverse took place in winter because low winter insolation cooled low latitudes more than high latitudes, which receive little winter insolation even today (Berger, 1978).A southern steepening of the summer LTG during the early-Holocene, however, represents a departure from the insolation trend (thin line, Figure 2f).Much of the change, prior to ca. 8 ka, corresponded with the declining influence of the Laurentide ice sheet on atmospheric waves over the Atlantic.The summer anomaly may indicate a thermodynamic response to the ice albedo anomaly (Morrill et al., 2018) rather than the year-round effect of the ice as a mechanical barrier to the flow of the jet stream (COHMAP, 1988).Paradoxically, the presence of the ice sheet thus combined with insolation anomalies to facilitate the northward expansion of oak at 10-8 ka.(Foster et al., 2006), (c) lake oxygen isotope records from the southern and northern regions (Grinnell Lake, New Jersey, Zhao et al., 2010;Fayetteville Green Lake, New York, Kirby et al., 2002), and (d) sea-surface temperature (SST) reconstructions from the Florida Strait (Lochte et al., 2020) and Labrador Shelf (Schmidt et al., 2012).

The Mid-Holocene Anomaly
The most significant departure from the insolation and ice sheet effects arose when a steep front developed in the north from 4.8 to 3.8 ka (Figures 2e,2f and 3e,3f).The front cooled the north but limited cold outbreaks to the south.It favored northern birch populations (Figure 2c), even as it enabled oak-hickory forests to expand into newly warmed central highlands (Figure 4a).The changes began as temperatures fell in the north after 6.2 ka (blue lines, Figures 2d and 3d) and were amplified when central temperatures reached their maximum at 4.8-3.8ka (black lines, Figures 2d and 3d).The millennial scale of this central temperature maximum has analogies at the other latitudes, such as the cluster of northern sites that warmed from 3.4 to 1.8 ka (Figures 2c and 3c).They indicate that more than one fluctuation in the LTG took place with a steep front returning to the north after 1.8 ka (Figures 2c-2f).
The engimatic absence of hemlock undoubtedly affected forest dynamics after 5.5 ka, but the similarity of the decline across latitudes (Figures 3b and 3c) contrasts with the north-south vegetation differences 700 years later (Figures 2b and 2c).Extensive droughts may have been important (Booth et al., 2012;Haas and McAndrews, 2000;Newby et al., 2014).They caused fluctuations in hemlock abundance at some sites (Marsicek et al., 2013) and are faithfully reconstructed from fossil pollen using our methods and data (Shuman et al., 2019).Other forest influences, such as wildfire and insect outbreaks, are not widely evident (Clark et al., 1996;Oswald et al., 2016Oswald et al., , 2020)), but local ecological interactions must have shaped landscape-scale changes (e.g., potentially explaining synchroneity of the decline with drought onset at one site but not another amid localized hemlock variability; Booth et al., 2012).Regardless, hemlock's demise did not drive all of the major forest changes.Abrupt changes in coastal oak and beech populations where hemlock is rare appear consistent with local cooling (Figure 4b; Foster et al., 2006).In the north, the simultaneous switch from hemlock to birch, another cold-tolerant tree (Figure 2c), raises a key question: why did oak replace hemlock and other cold-adapted trees like birch and beech in the central region (Figures 2b and 4a)-but not along the coast (Figure 4b) or in the north as it did in the early-Holocene (Figure 2c)?
The reconstructed LTG changes point to a plausible answer.Different trees systematically replaced hemlock by latitude because temperatures changed in opposite directions in a manner analogous to the NAO (Figures 1c-1e).Cooling forced the northern limit of abundant hemlock southward as central warming extended southern forests with little hemlock northward.Droughts amplified the pattern.Positive NAO phases similarly steepen the northern LTG today (Figure 1e), which is consistent with isotopic evidence of a steep northern front at 4.8-3.8ka (Figure 4c; Kirby et al., 2002).The NAO pattern extends warming from the mid-Atlantic region into the central study area despite little change in the southern study area (Figure 1d).Ocean cooling simultaneously develops from Labrador to Massachusetts (Figure S2 in Supporting Information S1), which would explain why a narrow coastal zone of oak-beech changes (Figure 4b) differs from the rest of the southern region (Figure 2a).
A contrast between low sea-surface temperatures (SSTs) from 4.5 to 3.5 ka on the Labrador Shelf (Lochte et al., 2020) and high SSTs in the Florida Strait from 4.7 to 3.3 ka (Schmidt et al., 2012) (Figure 4d) agrees with the patterns associated with a strong Atlantic pressure gradient today (Figure S2 in Supporting Information S1).However, a NAO reconstruction from Greenland indicates a negative phase at the same time (Olsen et al., 2012), which raises questions about whether the circulation changes involved have direct analogs to the NAO at monthly to annual scales.Regardless of the relationship, the mid-Holocene oscillation coincides with millennial anomalies in North Atlantic deep water flow, wind speeds, and atmospheric circulation (Giraudeau et al., 2010;Jackson et al., 2005;O'Brien et al., 1995).Evidence of unusual climates at the same time extends from paleosol development in Nebraska's Sand Hills (Miao et al., 2007) to isotopic anomalies in Africa (Thompson et al., 2002).The patterns may relate to intrinsic atmospheric variability, particularly if ocean-atmosphere or sea-ice interactions sustained the changes (Orme et al., 2021;Rigor et al., 2002;Thornalley et al., 2009).Alternatively, volcanic and solar variability may have been influential factors (Fletcher et al., 2013;Kobashi et al., 2017), but no clear linear correlation with external forcing is evident, potentially as expected (Renssen et al., 2006).

Conclusions
The Holocene temperature history of eastern North America includes changes in the steepness of the north-south temperature gradient, which generally responded in different directions to summer and winter insolation anomalies.Well-described Holocene forcing, however, does not explain an apparent millennial-scale latitudinal shift in

Figure 1 .
Figure 1.The fossil pollen sites (square symbols in a) span the gradient of modern mean annual temperatures in the northeast U.S. Additional maps represent modern variability in the latitudinal temperature gradient (LTG) in the study region (inset boxes in b-e).Maps of composite mean surface atmospheric pressures show years when the LTG was steepest (b) south of 42.75°N (1960, 1977, 1988, 1999, 2001) and (c) north of 42.75°N (1954, 1959, 1968, 1991, 2007).Lower maps represent (d) the associated composite mean temperature differences (LTG steep in north minus steep in south) and (e) the correlation coefficeient, r, between the North Atlantic Oscillation and mean annual temperature (T) in each U.S. climate division(Vose et al., 2014).Modern mean annual temperatures across the region (a) are based on 1991-2020 PRISM normals (PRISM, 2020); symbol colors indicate the geographical cluster containing each fossil site.See Supporting Information S1 for composite anomaly methodology.

Figure 2 .
Figure 2. Mean summer (June-August) temperatures as well as oak (Quercus) and birch (Betula) pollen percentages from (a) southern, (b) central, and (c) northern sites.Mean reconstructions for each cluster (bold lines) as well as the best resolved records within the 300-400 m elevation zone (thin lines) are shown in (d) alongwith the mean differences between clusters (e) and the north minus south difference in gradient steepness ("Gradient Diff." as the black line in (f).An inverted June insolation curve for 30°N(Berger, 1978) is normalized to the gradient difference scale for comparison (thin line, f).Vertical gray lines mark major time periods of change with the bold line at 5.5 ka indicating the median beginning of the regional hemlock decline (see Figure3).The vertical gray line at 4.8 ka is unlabeled.

Figure 3 .
Figure 3.As in Figure 2, but for mean winter (December-February) temperatures, hickory (Carya) and hemlock (Tsuga) pollen.A January insolation curve for 30°N (inverted, gray line) is normalized to the N-S gradient difference in f for comparison.

Figure 4 .
Figure 4. Representative oak, hemlock, and beech pollen records (a) from the central region (Guilder Pond and Little Pond, Royalston; Oswald et al., 2018) are compared with (b) pollen data from Deep Pond in coastal Massachusetts(Foster et al., 2006), (c) lake oxygen isotope records from the southern and northern regions (Grinnell Lake, New Jersey,Zhao et al., 2010; Fayetteville  Green Lake, New York, Kirby et al., 2002), and (d) sea-surface temperature (SST) reconstructions from the Florida Strait(Lochte et al., 2020) and Labrador Shelf(Schmidt et al., 2012).