Tectonics and climate are major contributors to the topographic evolution of mountain ranges. Here, we investigate temporal variations in exhumation due to the onset of Pleistocene glaciation in the Olympic Mountains (Washington State, USA). We present 29 new apatite and zircon (U-Th)/He ages (AHe and ZHe), showing a decrease in ages toward the interior of the mountain range for both thermochronometric systems. Young AHe ages (<2 Ma) can be found on the western side and the interior of the mountain range. Thermokinematic modeling of sample cooling ages suggests, that ZHe ages can be explained by an ellipse-shaped exhumation pattern with lowest/highest rates of ∼0.25 and 0.9 km/m.y. These rates are interpreted as tectonically driven rock uplift, where the pattern of rates is governed by the shape of the subducted plate. However, the youngest AHe ages require a 50–150% increase in exhumation rates in the past 2–3 m.y. This increase in rates is contemporaneous with Pliocene-Pleistocene alpine glaciation of the orogen, indicating that tectonic rock uplift is perturbed by glacial erosion.


The evolution of mountain topography (e.g., relief, mean elevation) is sensitive to variations in climate and tectonics that modulate the efficiency of various surface processes (e.g., Whipple, 2009). The onset of Pleistocene glaciation is hypothesized to have increased orogen exhumation rates, and significantly modified topography (e.g., Brocklehurst and Whipple, 2002; Ehlers et al., 2006; Valla et al., 2011; Glotzbach et al., 2013; Herman et al., 2013). Advances in low-temperature thermochronology and thermal modeling enable the quantification of spatial and temporal variations in exhumation (e.g., Braun, 2003). Here we test the hypothesis that enhanced Pleistocene glacial erosion can perturb the flux steady state of an orogen by increasing the erosional flux over million-year time scales.

We evaluate this hypothesis through an application to the tectonically active and glaciated Olympic Mountains located in Washington State, USA (Fig. 1A). This orogen is the exhumed portion of the Cascadia Subduction zone accretionary wedge (Tabor and Cady, 1978). Previous studies have suggested that exhumation rates have been largely constant since ca. 14 Ma, and that the orogen is in flux steady state, where accretionary and erosional fluxes are balanced (Brandon et al., 1998; Batt et al., 2001; Pazzaglia and Brandon, 2001). Largely unexplored in previous work is the potential transient effect of Pleistocene glaciation on the orogen-wide exhumation.

Here we complement previous work with new apatite and zircon (U-Th)/He ages (AHe and ZHe, respectively) from the Olympic Mountains (Fig. 1B) and compare them to predicted thermokinematic model ages to discriminate between different exhumation histories.


At the Cascadia subduction zone, the Juan de Fuca plate subducts beneath the North American plate and displays a three-dimensional (3-D) bend beneath the Olympic Mountains (Fig. DR1a in the GSA Data Repository1). The mountain range is dominated by two tectonostratigraphic units, which are separated by the Hurricane Ridge fault (HRF; Fig. 1; Fig. DR1b). The footwall represents the actual accretionary wedge, consisting of Eocene to Miocene marine sandstones and siltstones accreted to North America (Tabor and Cady, 1978; Brandon et al., 1998).

During the Pleistocene, alpine glaciers incised deep, wide valleys (e.g., Hoh, Queets, Quinault, and Elwha valleys; see Fig. 1A), locally forming piedmont-style glaciers approaching the Pacific Ocean (Thackray, 2001; Montgomery, 2002; Adams and Ehlers, 2017). Furthermore, the Cordilleran Ice Sheet surrounded the orogen in the north and east (Booth et al., 2003). Although the highest point (Mount Olympus, 2428 m) is located to the west of the range divide, most of the high topography is located east of the Elwha valley (Figs. 1A and 2). A strong gradient in modern precipitation exists, where the west side of the range receives 3–6 m/yr compared to the drier east (1–3 m/yr; see Fig. 2A and Fig. DR1c). Linked to the precipitation gradient, the Pleistocene equilibrium line altitude (ELA) increases from 1000 m on the west side of the divide to 1800 m on the east (Fig. DR1d).

Sedimentary rocks often contain thermochronometric age populations controlled by bedrock cooling histories. Burial and heating of sedimentary rocks (e.g., during metamorphism or subduction) reset thermochronometers, directly relating them to the exhumation of the collected bedrock sample. During exhumation, samples cool and pass through their closure temperatures: ∼240 °C for zircon fission-track (ZFT), ∼180 °C for ZHe, ∼100–120 °C for apatite fission-track (AFT), and ∼60–70 °C for AHe (Reiners and Brandon, 2006). Published thermochronometer ages for the Olympic Mountains are compiled in Figure DR3. Young AFT ages (2.3–3.9 Ma) are present on the western side of the mountain range, but the only four reset ZFT samples (13–14 Ma) are located to the east of Mount Olympus. Youngest AHe ages (2.0–2.5 Ma) can be found in the vicinity of Mount Olympus. High exhumation rates of 1 km/m.y. are suggested on the west side of the range (Fig. DR1e), which then decrease toward zero near the coast of the peninsula (Brandon et al., 1998; Pazzaglia and Brandon, 2001).


We collected 30 bedrock samples for (U-Th)/He thermochronometric dating at an elevation of ∼400 m (Fig. 1B; Fig. DR2). This equal-elevation sampling approach allows for direct comparison of ages, without the need to correct ages for different sample elevation. Twenty-nine (29) of these samples yielded datable apatite crystals (e.g., suitable crystal shape, no inclusions) and 27 samples were selected for ZHe dating. A map of AHe and ZHe ages (Fig. DR3), a description of analytical methods, and data tables are provided in the Data Repository.

Reset AHe ages (27 samples) range from 15.0 to 1.5 Ma, but two samples at the east coast contain un-reset ages (Fig. DR3). Eleven samples are younger than or equal to the onset of glaciation (ca. 2–3 Ma), which are all located to the west of the Elwha valley. Only 14 samples have reset ZHe ages (4.8–14.3 Ma), and they are all located in the deeply incised valleys within the high-topography part of the range. The youngest ZHe ages (4.8 Ma and 5.9 Ma for samples OP1533 and OP1515, respectively) can be found in the headwaters of the Hoh and Elwha valleys.

Swath profiles of our data perpendicular (A-A′) and parallel (B-B′) to the range divide are shown in Figure 2 (for location, see Fig. 1B). In both swath profiles, AHe and ZHe ages decrease or change from un-reset to reset toward the center of the range. For swath B-B′, the high-topography part of the range overlaps with the area of young cooling ages (kilometers 20 and 70 of the swath). However, in swath A-A′, the area of reset ZHe samples and young AHe ages (<2.5 Ma) is offset from the highest topography and shifted toward the west.


Converting thermochronometric ages into exhumation histories requires estimation of the geothermal gradient over time. For this, we use the 3-D thermokinematic model Pecube (Braun, 2003), integrated with previous ages and our new thermochronometer ages. Although the importance of frontal accretion and horizontal velocities in the Olympic Mountains has been explored by previous studies (Batt et al., 2001), we only consider vertical movement in our model simulations. Further details and discussion about the justification of this approach, as well as explanations about the modeling, can be found in the Data Repository.

Our thermokinematic model uses the present-day topography of the Olympic Mountains as input. Following Brandon et al. (1998), we start exhumation and development of topography at 18 Ma, and reach a steady state at 14 Ma. Orogens situated in an orogenic syntaxis (like the Olympic Mountains) are predicted to show a “bull’s eye” pattern of exhumation (Bendick and Ehlers, 2014). Therefore, we choose an ellipse-shaped exhumation pattern for our model simulations (Fig. 1B), matching mean elevation and relief (Fig. DR7), which are suggested to correlate with rock uplift (Adams and Ehlers, 2017). The ellipse pattern is defined by a minimum rate at the edge of the ellipse and a maximum rate at the centroid of the ellipse. Our first objective is to find a spatial pattern of constant, long-term exhumation rates. To achieve this, we vary the location and size of the ellipse, and the gradient in exhumation rates, to minimize the misfit between modeled and observed ages (from literature and our data) with a reduced χ2-test (Fig. DR6). Our second objective is to test the hypothesis that the exhumation rates increased in Pliocene-Pleistocene times due to enhanced glacial erosion. To thoroughly investigate plausible increase scenarios, we take the preferred ellipse from the previous step and increase the rate of the entire model domain each at six different time steps by seven different magnitudes. From the resulting 42 combinations, we find the best-fit time and magnitude by comparing our and published ages with the modeled ages using a reduced χ2-test (results are shown in Figure DR8 and Table DR5). Finally, for comparison with our new data, modeled ages are extracted along the swaths A-A′ and B-B′ (Fig. 1B) from elevations between 200 and 500 m and displayed in Figure 2.

Our preferred ellipse for temporally steady exhumation rates has diameters of 70–60 km and the centroid of the ellipse is located in the core of the mountain range. Rates of 0.25 km/m.y. are found at the edges and increase to 0.9 km/m.y. in the center (Fig. 1B). Using this ellipse and keeping the rates constant throughout the entire model duration results in ages (black solid lines in Figs. 2B, 2C, 2E, and 2F) reproducing the general age pattern for AHe and ZHe. Abrupt changes in modeled ages (particularly for ZHe) are caused by the strong gradient in exhumation rates. However, for constant-rate models, the modeled AHe ages are anomalously old compared to observed ages, indicating that an increase in rates is required. We find a wide variety of transient exhumation histories, where rate increases occurred between 10 and 2 Ma. These transient exhumation histories produce improved statistical fits to the data (total χ2 values of 50–54; Fig. DR8) than the constant-rate simulations (lowest total χ2 value of 66; Fig. DR6). A timeline summarizing possible exhumation histories from our transient simulations is shown in Figure 3.


Our new data and modeling allow us to resolve the long-term spatial pattern of exhumation in the Olympic Mountains and to identify an increase in exhumation that occurred around the same time as the onset of Pliocene-Pleistocene glaciation. The strong spatial gradient in thermochronometer ages (e.g., Fig. 2) requires a proportional spatial gradient in exhumation, and our best-fit model suggests an elliptical exhumation pattern with minimum and maximum rates of 0.25 and 0.9 km/m.y. While the range of exhumation rates agrees with previous findings (Brandon et al., 1998; Pazzaglia and Brandon, 2001), our spatial pattern of exhumation rates is different. The highest exhumation rates occur in the high, rugged core of the range, encompassing the headwaters of Hoh, Queets, Quinault, and Elwha rivers (Fig. 1B), as also suggested by Adams and Ehlers (2017) based on topographic analysis. We interpret that the general ellipse-shaped exhumation pattern is imposed by the geometry of the subducted plate, similar to other syntaxial orogens with focused exhumation (Bendick and Ehlers, 2014; Falkowski and Enkelmann, 2016; Lang et al., 2016). The required bending of the subducted plate and slab convexity is well known in the Olympic Mountains (e.g., Crosson and Owens, 1987; see also Fig. DR1a).

A temporally constant exhumation rate may explain the observed cooling ages of higher-closure-temperature systems (AFT, ZHe, ZFT), but transient exhumation histories provide much better statistical fits for ZHe, and for the AHe ages. The χ2-test also reveals that there is no unique exhumation history solution, and several transient scenarios produce equally good fits (Fig. DR8). Generally speaking, the older the increase in rates is, the smaller is the required amount of increase (e.g., 30–50% rate increase at 7 Ma compared to 100–150% at 2 Ma). In the absence of other tectonic or climate mechanisms that might increase exhumation rates at older times (e.g., >5 Ma), increased erosion due to glaciation during the Pliocene-Pleistocene is the most plausible mechanism for explaining the young increase scenarios at 2–3 Ma with an increase in rates by 50–150%. These transient scenarios (colored envelopes in Figs. 2B and 2E) result in AHe ages that are 2–4 Ma younger compared to constant-rate predictions, thereby improving the fit to observed ages. However, temporal constraints for the onset of alpine glaciation in the Olympic Mountains are ambiguous, and areas nearby experienced glaciations earlier than 2–3 Ma (Fig. 3). Indeed, an earlier increase in exhumation at 4 or 5 Ma (but with a smaller magnitude) is permissible (colored dashed lines in Figs. 2B, 2C, 2E, and 2F). The ELA corresponds to the area of most effective glacial erosion (e.g., Montgomery, 2002) and is lowest on the western side of the range (Fig. DR1d). Thus, the effect of glacial erosion on exhumation is expected to be strongest on the west side of the mountain range. Here, the AHe ages from three samples (OP1527, OP1529, and OP1532) support this hypothesis and locally suggest an increase in rates by 150–200% at 2–3 Ma (gray zones in Fig. 2B). This area could also have experienced a glacial-related change in topography, causing younger AHe ages (e.g., Ehlers et al., 2006). Other reasons that might explain a mismatch between modeled and observed ages include a deviation from the assumed, perfect ellipse-shaped exhumation pattern or an even more complex transient exhumation history.

In conclusion, the Olympic Mountains are the product of both spatial and temporal variations in exhumation rates. While the spatial pattern of exhumation is governed by the tectonic setting, the temporal variation was caused by Pliocene-Pleistocene glaciation. The magnitude of exhumation rate increase is similar to other neighboring orogens influenced by Pleistocene glaciation (e.g., Ehlers et al., 2006). Given our observed temporal increase in rates, and the increase in material flux and the change in deformational style of the wedge during the Quaternary (Adam et al., 2004), the proposed flux steady state of the mountain range (Batt et al., 2001) is questionable.


This work was supported by a European Research Council (ERC) Consolidator Grant (615703) to Ehlers. We thank Holger Sprengel, Bill Baccus, Jerry Freilich, Roger Hoffman, and the Olympic National Park rangers for their help and assistance during field work. We thank Mark Quigley for editorial handling of the manuscript, and Doug Burbank, Jim Spotila, and an anonymous reviewer for constructive reviews.

1GSA Data Repository item 2018161, supplementary figures, details about methods and modeling, and analytical results (in Excel tables), is available online at http://www.geosociety.org/datarepository/2018/ or on request from editing@geosociety.org.
Gold Open Access: This paper is published under the terms of the CC-BY license.