Holocene relative sea-level changes and glacial isostatic adjustment of the U . S . Atlantic coast

The first quality-controlled Holocene sea-level database for the U.S. Atlantic coast has been constructed from 686 sea-level indicators. The database documents a decreasing rate of relative sea-level (RSL) rise through time with no evidence of sea level being above present in the middle to late Holocene. The highest rates of RSL rise are found in the mid-Atlantic region. We employ the database to constrain an ensemble of glacial isostatic adjustment models using two ice (ICE-5G, ICE-6G [global ice sheet reconstructions]) and two mantle viscosity (models VM5a,VM5b [VM—radial variation of viscosity in the sublithospheric mantle]) variations to assess whether the spherically symmetric viscoelastic models are able to survive intercomparison with a more refined database of postglacial RSL history. We identify significant misfits between observations and predictions using ICE-5G with the VM5a viscosity profile. ICE-6G provides some improvement for the northern Atlantic region, but misfits remain elsewhere. Decreasing the upper mantle and transition zone viscosity by a factor of 2 to 0.25 × 10 21 Pa s (VM5b) removes significant discrepancies between observations and predictions along the mid-Atlantic coastline, although misfits remain in the southern Atlantic region. These may be an indication of the importance of laterally heterogeneous viscosity in the upper mantle.


INTRODUCTION
Relative sea-level (RSL) change is related to the redistribution of mass from ice sheet growth and decay, inducing isostatic compensation of the underlying solid Earth.Our understanding of current rates of sea-level rise from tide gauge (e.g., Church and White, 2006) and satellite (e.g., Cazenave et al., 2009) data, and of the ongoing mass loss from the major ice sheets by the Gravity Recovery and Climate Experiment (GRACE) (e.g., Velicogna and Wahr, 2006) requires correction for glacial isostatic adjustment (GIA) effects that are both calibrated to, and independently tested by, observations of former sea levels.Holocene RSL data are used to infer mantle viscosity (e.g., Mitrovica and Peltier, 1995) and lithospheric thickness (e.g., Tushingham and Peltier, 1992).The Atlantic coast of the United States is a key region for the comparison of model predictions and sea-level observations because it provides an independent constraint for GIA models such as ICE-5G VM5a (Peltier and Drummond, 2008) that are tuned to data sets from Canada, Fennoscandia, and Barbados (Peltier and Fairbanks, 2006).It is signifi cant that the GIA of the U.S. Atlantic coast is dominated by the collapse of the largeamplitude proglacial forebulge of the massive Laurentide ice sheet.
The earliest GIA models (Clark et al., 1978) did not fi t the observational data from the U.S. Atlantic coast when the viscosity of the mantle was assumed to be independent of depth (Cathles, 1975).While the VM1 (VM is radial variation of the viscosity of the sublithospheric mantle) model improved the fi t between observations and model predictions (Tushingham and Peltier, 1992), a better agreement was achieved with the combination of ICE-4G (global ice sheet reconstructions) and the more complex VM2 viscosity profi le (Peltier, 1996).Since the publication of Peltier (1996), however, there have been advances in the reconstruction of RSL (e.g., Horton et al., 2009) and the development of improved ice models, including ICE-5G (Peltier, 2007) and ICE-6G (Peltier, 2010).The incorporation of rotational feedback in GIA models (e.g., Wu and Peltier, 1984) is an especially important recent advance.Here we focus on whether a new sea-level database is able to confi rm the good quality of fi t between observations and predictions previously obtained using simpler models, and reveal signifi cant systematic deviations.

CONSTRUCTION OF A SEA-LEVEL DATABASE
The sea-level database is created from published and unpublished samples of organic sediment (salt and fresh-water marshes) and shells of marine gastropods and foraminifera.The samples are converted into sea-level index points (see Supplementary Methods in the GSA Data Repository 1 ) when they meet three criteria: (1) the location of the sample is established to within 1 km; (2) the age of the sample is calibrated to sidereal years using the latest 14 C calibration curve; and (3) a relationship between the sample and a known tidal level (i.e., indicative meaning; van de Plassche, 1986) is identifi ed.For samples that cannot be directly related to former sea level, we produce marine (e.g., marine shells) and terrestrial (e.g., freshwater peat) limiting data from samples that must have been deposited below and above mean sea level, respectively.Every sample has an error calculated from a variety of factors inherent to sealevel research (van de Plassche, 1986).We minimized the infl uence of compaction of sediment by excluding intercalated index points (organic sediments that were underlain and overlain by different sedimentary units) and subdivide the remaining into "base of basal" and "basal" (e.g., Horton and Shennan, 2009).To account for spatial variations, the database was subdivided into 16 geographical regions (1-16 in Fig. 1) based on a combination of the availability of data, the distance from the former center of the Laurentide Ice Sheet, and GIA reconstructions.The addition of new RSL data may allow further subdivision.For all observational data (radiocarbon dates, calibrated age ranges, RSL reconstructions with associated errors, and references to the original publications), see the Data Repository.

PREDICTIONS OF RSL
The model analyses are based on the full gravitationally self-consistent form of the GIA theory and include the effects of rotational feedback (e.g., Milne and Mitrovica, 1996;Peltier, 1998;Wu and Peltier, 1984).The RSL predictions are based on the ICE-5G (Peltier and Drummond, 2008) and ICE-6G (Peltier, 2010) ice models.The ice models are coupled to the VM5a viscosity model (Peltier and Drummond, 2008) that reduces the misfi t between predicted and observed horizontal motions of the North American plate (Argus and Peltier, 2010).VM5a includes a 100-km-thick lithosphere, consisting of a 60-km-thick elastic upper layer, beneath which there is a 40-km-thick layer with a viscosity of 10 22 Pa s.Here we revise the Earth model by reducing the viscosity of the upper mantle and transition zone above the 660 km phase transformation from 0.5 × 10 21 Pa s (VM5a) to 0.25 × 10 21 Pa s (VM5b) (Fig. DR2 in the Data Repository).This is similar to the value calculated by Wolf et al. (2006) of 0.32 × 10 21 Pa s but represents a 50% reduction compared to the viscosity of 0.53 × 10 21 Pa s provided by Paulson et al. (2007).However, Paulson et al. (2007) noted that a lower upper mantle viscosity could be supported if associated with a stronger lower mantle viscosity.In Tushingham and Peltier (1992) it was shown that reducing the viscosity of the upper mantle and transition zone substantially increases the width of the proglacial forebulge, and therefore raises the predicted RSLs along the U.S. Atlantic coast without increasing the number of free parameters.VM5b continues to fi t RSL data from northern Canada because they are insensitive to upper mantle and transition zone viscoelastic structure (Peltier, 1998).

HOLOCENE SEA-LEVEL OBSERVATIONS
The new Holocene RSL database for the U.S. Atlantic coast consists of 342 sea-level index points, 189 marine limiting dates, and 155 terrestrial limiting dates (Fig. 2A).The RSL database previously used to constrain GIA models (Peltier, 1996) contained fewer index points (n = 175) and marine limiting data (n = 85), but a greater number of terrestrial limiting data (n = 395).The increase in the number of index points in the new database is due to both the addition of new data (e.g., Horton et al., 2009;Miller et al., 2009) and the reinterpretation of terrestrial limiting dates as index points on the basis of the macrofossil and microfossil sea-level indicators.The elevation errors are index-point specifi c, in contrast to the standard vertical error term employed in the previous database.The new database has good temporal coverage from 6 ka to present; however, only 7% of the index points are older.The early Holocene record is primarily defi ned by limiting data.Spatially, RSL is well constrained by index points between Maine and South Carolina, although there is an absence of index points in Georgia and the Atlantic coast of Florida.The validation of observations strongly affects the interpretation of a GIA model.For example, the New Jersey predictions must now plot through index points, rather than being below a host of samples previously interpreted as limiting data (Fig. 2B).The North Carolina site demonstrates the importance of the addition of new data (Horton et al., 2009) to constrain models of GIA (Fig. 2C).
Analysis of the full database from eastern Maine to southern South Carolina (regions 1-16 in Fig. 3) demonstrates that RSL has not risen above present during the middle and late Holocene.In Maine and northern Massachusetts, limiting data suggest that RSL may have dropped from above present prior to a slowstand between 11.5 and 7.5 ka (e.g., Kelley et al., 2010).Rates of RSL change were highest during the early Holocene and have been decreasing over time, due to the exponential form of the GIA process following deglaciation and the reduction of ice equivalent meltwater input from 7 ka onward.The maximum rate of Holocene RSL rise (e.g., ~15 m since 6 ka) occurred in New Jersey and Delaware (regions 8−10), coincident with the area of greatest ongoing GIA related subsidence.The RSL histories of the northeastern Atlantic region (regions 1−5) are constrained by sea-level index points from 7 ka to present; all fi ve areas reveal a rise in sea level of <10 m since 6 ka.Index points from southern North Carolina to southern South Carolina (regions 14−16) similarly support a rise in RSL of <10 m during the past 6 k.y.

IMPLICATIONS FOR GLACIAL ISOSTATIC ADJUSTMENT MODELS
The ICE-5G VM5a model is in good agreement with the Holocene RSL observations in eastern (region 1) and southern Maine (region 2) for the past 6 k.y.The model is above marine limiting dates from southern Maine between 8 and 11 ka.For the remaining study areas (regions 3−16), the model fi ts the observations in the late Holocene (0-3 ka), but with increasing age there is a systematic disagreement between the model and the data.The misfi t is most pronounced along the mid-and southern Atlantic coastlines (New York, 6, to southern South Carolina, 16) with observations of RSL ~10 m higher than model predictions at 6 ka.The predictions are invalidated by marine limiting dates at southern Massachusetts (4), New Jersey (8), Chesapeake Bay (11), and northern North Carolina (13).ICE-6G VM5a is an improvement over ICE-5G VM5a for northern Massachusetts to New York (regions 3−6).However, the model now underpredicts RSL in eastern Maine (1) and overpredicts in southern Maine (2).There is little difference between ICE-5G and ICE-6G from Long Island (7) to southern South Carolina (16), because both ice models have similar mass and cover exactly the same surface area of the North American continent.
For both ice models, decreasing the upper mantle viscosity in VM5a to produce VM5b results in a considerable improvement (Fig. DR3) in the quality of fi t along the U.S.  Atlantic coast in the area of greatest GIA-related subsidence (regions 4−12).Earlier model iterations that also increased the upper mantle/lower mantle contrast ratio from 1:1 to 1:4 resulted in a similar decrease in variance (Tushingham and Peltier, 1992).However, by reducing the value of both the upper mantle and the transition zone viscosity, we may have signifi cantly violated the model fi t to the McConnell (1968) spectrum of Fennoscandian rebound.This would suggest that we have directly detected a requirement for lateral viscosity variation in the upper mantle from two independent data sets associated with the near-fi eld Fennoscandia region and intermediate-fi eld U.S. Atlantic coast (e.g., Kaufmann and Lambeck, 2002).It has previously been noted that one-dimensional viscosity structure derived from global postglacial rebound observations is heavily weighted toward the viscosity structure beneath the loaded region and not the global average (Paulson et al., 2005).
While the VM5b viscosity profi le results in a signifi cant improvement, there remain two outstanding issues.First, the incorporation of the VM5b viscosity profi le with ICE-5G predicts late Holocene highstands of RSL in the northeastern United States (regions 1−3), because the softening of the upper mantle and transition zone produces a time-dependent shift of the boundary between uplift and subsidence.The highstand is removed at eastern Maine (region 1) by using ICE-6G because of a change in the thickness of the proximal ice load.Therefore, the remaining highstands in southern Maine (region 2) and northern Massachusetts (region 3) may also be eliminated through further thickening of the ice load in proximity to these two locations.Alternatively, a slight increase in the thickness of the elastic lithosphere would accomplish the same improvement of fi t (Tushingham and Peltier, 1992).
Second, and more important, are the continuing misfi ts between the VM5a/VM5b RSL predictions and observations in the southern Atlantic region (regions 13-16).We must consider that changes to other aspects of the Earth model may be necessary to fi t the data, including lithospheric thickness (e.g., Tushingham and Peltier, 1992) and/or the incorporation of lateral heterogeneity in the mantle (e.g., Wu, 2006) and lithosphere (e.g., Wang and Wu, 2006).In Tushingham and Peltier (1992), it was demonstrated that changes in lithospheric thickness from 71 to 245 km have little effect on the variance between model predictions and sea-level observations during the Holocene along the U.S. Atlantic coast.Three-dimensional mantle viscosity may be required to incorporate the effects of the subduction of the Farallon plate 80 m.y.ago, which has now penetrated the top half of the lower mantle beneath the Atlantic (e.g., Bunge and Grand, 2000).Changes to the viscosity of the upper part of the lower mantle in the one-dimensional model (e.g., Davis and Mitrovica, 1996) are usually ruled out by the strong constraints of the RSL data from Canada.However, Wang and Wu (2006) produced a model supported by seismic velocity anomalies from global tomography models that showed lateral variations at all depths of the mantle performed better than a model with no lateral variations.
RSL observations may also be subject to changes from processes such as compaction, tidal range variations, and tectonics (e.g., Shennan and Horton, 2002).The database shows no discernable difference between basal and base of basal samples; compaction would only lower the elevation of index points, artifi cially improving the fi t to the model.If the tidal range has not remained constant through time, sea-level chronologies based upon tidelevel indicators will differ from true sea level.Modern tidal range along much of the southeastern Atlantic coast is <2.5 m.Therefore, an unrealistic tenfold increase in tidal range is needed to lower the elevation of sea-level observations and remove the misfi t at 7 ka.The rate of tectonic uplift necessary to achieve a fi t between model and data (>1 mm yr -1 ) disagrees with global positioning system (Sella et al., 2007) and tide-gauge (e.g., Douglas, 1991) observations.Uplift rates of this magnitude are more typically associated with seismic coastlines such as Cascadia (e.g., Burgette et al., 2009).There is evidence for tectonic activity in North Carolina and South Carolina.However, the calculated long-term uplift rate of ~0.02-0.05mm yr -1 along the Orangeburg Scarp (South Carolina) (Dowsett and Cronin, 1990) is not suffi cient.The rate of uplift along  the Cape Fear Arch (North Carolina) is poorly constrained (0.14-1.8 mm yr -1 ) and may be highly localized (Marple and Talwani, 2004).Rates of change due to mantle upwelling and downwelling show no spatial pattern along the U.S. Atlantic coast and are likely to be <0.1 mm yr -1 (Conrad and Husson, 2009).It is unfortunate that there are currently no RSL data available from the Atlantic coast of Florida or Georgia to allow us to ascertain whether the misfi t between observations and predictions is a regional phenomenon.This limitation currently prevents us from confi rming that the internal viscoelastic structure must be laterally heterogeneous on a scale that infl uences ongoing GIA.

ACKNOWLEDGMENTS
Figure 1.Location map of United States Atlantic coast showing 16 study areas from Maine to South Carolina.Numbers correspond to reconstructions in Figure 3. ME-Maine; MA-Massachusetts; RI-Rhode Island; CT-Connecticut; NY-New York; NJ-New Jersey; PA-Pennsylvania; DE-Delaware; MD-Maryland; VA-Virginia; NC-North Carolina; SC-South Carolina; GA-Georgia; FL-Florida.

Figure 2 .
Figure 2. Age-altitude plots of relative sealevel (RSL) observations for Peltier (1996) and new database.A: Highlighting all index points; inset is histogram showing temporal distribution.B, C: Differences in New Jersey and northern North Carolina.Index points in new database are plotted as boxes with 2σ vertical and calibrated age errors, whereas in previous database only age errors are 2σ.