Post-glacial Rebound in the Okanogan Valley?

The Okanogan Valley runs 100km due north from Brewster, WA to the Canadian border; it provides a near perfect longitudinal profile through the Okanogan Lobe of the Cordilleran Ice Sheet. The "kame terrace" system preserved there is extensive and accessible. I call them kame terraces for simplicity; their origin is certainly more complicated. Kame terraces are accumulations of glaciofluvial sediment, glaciolacustrine sediment, and till along the margin of a glacier, between the bedrock sides of the valley and the slug of ice the (mostly) fills it. Their tops tilt gently down-valley governed by local base level determined by the ice body itself. The deposits are the result of ice-marginal sedimentation, which is a bit different from what happens in moraines or on outwash plains.

Recently published articles proposing isostatic rebound (surface uplift) following deglaciation in the Upper Grand Coulee area (Pico et al., 2022; Waitt et al., 2021) have raised some questions in my mind:

The Okanogan Valley terrace train, in part glaciolacustrine, stretches north for 100km and is incredibly well preserved. Post-glacial crustal rebound can be up to ~1/3 the thickness of the ice. If several meters of rebound have occurred at the glacial terminus since the retreat of ice, as is purported for Upper Grand Coulee, then as you go north (thicker ice, greater weight), the uplift signal should strengthen (i.e., like for the Puget Lobe). Elevation should increase and slopes should steepen with rebound.

- Is there tilt on the glacial lakebeds in Okanogan Valley or similar deposits in nearby valleys?

- Is the valley long profile oversteepened?

- Does the valley's glacial terrace system (relict surfaces) show more tilting on its north end or does the profile resemble the typical concave "equilibrium" profile of the modern Okanogan River?

- The topic of tectonic signals in river long profiles has received substantial attention in the literature, especially in Geomorphology and JGR Earth Surface (i.e., Whipple, Merritts, Snyder, Demoulin, Sklar and Dietrich, Seidl, Burbank and Anderson, Church, Whittaker, Montgomery, etc.).

Gridded ice thickness values in meters for the Cordilleran Ice Sheet at its farther-advanced position, ~14 ka (James et al., 2000, Fig. 4). Thickness values for the Okanogan Valley area range between 60 and 980, so 10 to 100 m of post-glacial rebound is a reasonable expectation if the 1/3 ice thickness rebound rule is used. Ten meters of rebound is in theory a large enough uplift signal to be seen in long profiles of the valley's glacial terraces. Circles reference ICE-3G modeling by Tushingham and Peltier (1991).

An analysis of GPS crustal velocity data for North America defined a northern region exhibiting significant post-glacial isostatic rebound (red arrows) and a southern region showing no rebound or subsidence (yellow arrows). The two regions are separated by a zero velocity "hingeline" (green line), which crosses most of the continent, but stops short of the Cordilleran Ice Sheet at the Rocky Mountain front. Note where the Okanogan Valley falls on the map (red circle) and how modest the modeled rebound signal is (blue arrows) in comparison to nearly everywhere north of the hingeline. Post-glacial isostatic adjustment in the Pacfic Northwest seems to operate differently than areas to the east. Perhaps ice-loading histories are different or mantle viscosity outboard of the craton differs from that below it. Map modified from Sella et al. (2007, Fig. 1). IGb00 is a global GPS elevation datum.

Vertical crustal motion modeled from gravity measurements (Paulson et al., 2007). Blue and purple areas indicate uplift (rise) following removal of Pleistocene ice sheets. Yellow and red areas indicate subsidence (fall). The Okanogan Valley falls within an area of modest uplift near the zero velocity contour.

Black contours are predicted uplift contours from Sella et al. (2007). While the observation data (red and blue arrows) don't extend far enough west, the Okanogan Valley would fall just north of the predicted zero velocity contour. Modified from Van der Wal et al. (2009, Fig. 4).

A long profile 1km in length shows terraces on one side of the Okanogan Valley near Omak, WA. Terrace flats were defined as "nearly flat" pixels (slope < 13%) using ArcGIS Spatial Analyst and a 4km-wide swath of 10m DEM data centered on the river centerline. Pixels were converted to x,y points and exported as .dbf to facilitate easy plotting in a spreadsheet. Data in the plot shows all pixels that met our "nearly flat" criteria up to 2km away from the river (west bank, river right). Relative elevation (y) is plotted against downstream distance (x). Color ramp corresponds with distance away from the river (red = close, yellow/tan = far). Smaller symbols would be preferable, but we didn't figure that out.

A DEM-based profile of a flight of glaciofluvial and kame terraces along the Okanogan Valley near Omak, WA. A number of terrace levels can be seen. Recall this is only a 1km-long profile. Artifacts such as flat spots in stream channels descending distant hillslopes are flagged.

My new work. Extracting terrace flats from lidar data covering the Okanogan Valley. The image above is 1-km long segment of the valley profile. X-axis is elevation. Y-axis is downstream distance. Dots represent pixels with nearly horizontal slopes. Only one side of the valley is shown. We can see a low outwash terrace system. A discontinuous mid-level kame terrace is also clear. Remnants of higher terrace surfaces and low-angle alluvial fans that grade to them are present, too. Diagonal features are gullies (flat spots along gullies are common).

A few years ago, when I was teaching for Boise State Geosciences Department, I had my students delineate kame terrace flats from 10m DEM data. See below for a description of that work. Today, high resolution lidar data area available for the Okanogan Valley, which would make terrace mapping much more precise.

Terraces from a DEM - Each data point is a pixel value created from a simple calculation for 'relative elevation' made on the 30m DEM. Thus, the values are not actual elevation (but scale linearly with it) and many of the pixels in the original DEM were filtered out (pixels with slopes greater than about 10%, which is very flat). Each point shown in the chart satisfies the filter criteria (nearly level slope, within 2km from river, river right only). The color ramp represents distance away from the river (red=close, yellow = distant). My annotations/interpretation are in pen.

The blurriness in some terraces in the above example is the result of four factors other than the DEM resolution (10m x 10m pixels):

1.) The size of the symbol we used to plot the data could only be made so small (we wanted to make them smaller). We used open circles here, but horizontal dashes might be better. We could vertically stretch the display to change how things look a bit, too.

2.) The terraces slope gently down-valley (to the left) as well as toward the river (in a plane out of the screen).

3.) The surfaces of glacial terraces like these often undulate a bit (kettles, drainages, etc).

4.) Terraces formed primarily by deposition of sediment by flowing water in the presence of ice. Such deposition is complex with local cuts and fills, which express as variations in the y-axis (elevation).

Since we know each charted point represents a flat-lying pixel, we can interpret individual terrace surfaces in a coarse way (lump all the blue pixels in the example) or split them into individual terraces (black lines drawn on blue pixels).

This analysis is probably best done for long river valleys using a combination of ArcGIS (DEM clipping, calculations) and R (looping, charting). There are other technologies that would work, too. Terrace interpretation is a human task no matter how you do it.

This method could readily be used to complete a regional analysis of the terrace system along the former margin of the Cordilleran Ice Sheet in Washington, including the so-called "Great Terrace".

Kiver and Stradling partially tackled the kame terrace problem back in the 1990s using 24k topo maps.

The Great Terrace near Chelan, WA from I.C. Russell's USGS reconnaissance report of 1893. How does it relate to extensive terrace trains in the valley to the north?

A rapid topographic profile I sketched in Google Earth follows a high level surface that continues from the Columbia Valley to the Okanogan Valley. The profile (dashed blue line) seems to show a match between the highest terraces in the lower Okanogan Valley and the the Great Terrace of the Columbia at its northern end. Not sure. More to come.

References

James, T.S.; Clague, J.J.; Wang, K.; Hutchinson, I., 2000, Postglacial rebound at the northern Cascadia subduction zone, Quaternary Science Reviews, v. 19, p. 1527-1541

Russell, I.S., 1893, A geological reconnoissance in central Washington, USGS Bulletin 108

Russell, I.S., 1898, The Great Terrace of the Columbia and other topographic features in the neighborhood of Lake Chelan, WA, American Geologist, v. 22

Paulson, A.; Zhong, S.; Wahr, J., 2007, Inference of mantle viscosity from GRACE and relative sea level data, Geophysical Journal International, v. 171, p. 497–508

Pico, T.; Lamb, M.; Larsen, I.; David, S.; Mix, A., 2022, Glacial isostatic adjustment directed incision of the Channeled Scabland by Ice Age megafloods, Proceedings of the National Academy of Sciences, v. 119

Tushingham, A.M.; Peltier, W.R., 1991, ICE-3G: a new global model of late Pleistocene deglaciation based upon geophysical predictions of post-glacial relative sea level change, Journal of Geophysical Research, v. 69, p. 4497-4523

Waitt, R.B., Atwater, B.F., Lehnigk, L., Larsen, I.J., Bjornstad, B.J., Hanson, M.A., O'Connor, J.E., 2021, Upper Grand Coulee: New views of a channeled scabland megafloods enigma. enigma, in Booth, A.M., and Grunder, A.L. (editors), From Terranes to Terrains: Geologic Field Guides on the Construction and Destruction of the Pacific Northwest: GSA Field Guide 62, p. 245–300

Van der Wal, W.; Braun, A.; Wu, P.; Sideris, M.G., 2009, Prediction of decadal slope changes in Canada by glacial isostatic adjustment modeling, Canadian Journal of Earth Sciences, v. 46, p.. 587-595