An Introduction to GIS-Based Techniques for
Reconstruction of Quaternary Landscapes
Digital Reconstructions of Quaternary Landscapes
During glacial periods, the formation and expansion of ice sheets results in loading of the Earth’s crust, causing terrain subsidence at ice sheets and in broad regions adjacent to ice margins (e.g., Walcott, 1970; Peltier, 1974; Peltier and Andrews, 1976; Peltier, 1985, 2004). Later reductions of ice mass during deglaciation cause terrain emergence through glacio-isostatic rebound (e.g., Walcott, 1972; Farrell and Clark, 1976; Clark et al., 1978; Andrews, 1989; Dyke and Peltier, 2000). The extent to which the surface of a region of interest has been uplifted since a given time can be estimated if preserved strandline features, formed at a time of interest and in association with an extensive body of water, are widely distributed across the region. The spatial pattern of emergence since a given point in time can be approximated by plotting the modern elevations of the once-horizontal but now differentially raised strandlines of a particular age on a map. If sufficient strandline data are available, the modern geometry of a former marine or lacustrine surface can be approximately defined using isobases, which are curves that demarcate zones of equal glacio-isostatic rebound (e.g., Johnston, 1946; Andrews, 1966; McCann and Chorley, 1967; Smith et al., 1969; Andrews, 1970; Blake, 1970; Andrews et al., 1971; Dyke, 1974). Isobases are contours that describe the relative morphometry of time-equivalent crustal deformation resulting from glacial unloading (Andrews, 1970), although the magnitudes of isobase values can also be a function of other factors including changes in sea level (e.g., Walcott, 1972; Andrews, 1989).
Computer-based reconstruction of a glacial or post-glacial landscape allows for the generation of paleo-topographic databases and for quantitative analysis of the nature of ancient water bodies and drainage networks, permitting the establishment of numerical constraints on environmental parameters of relevance to the study of geomorphology, oceanography, and climatology. Such reconstructions can be developed for a particular time of interest by subtracting interpolated isobase values from modern elevations. The resultant paleo-topographic databases can be used to characterize past environments and investigate past geological, hydrological, and climatic processes (e.g., Gareau et al., 1998; Lewis and Gareau, 2001; Leverington, 2002; Leverington et al., 2000, 2002ab; Schaetzl et al., 2002; Leverington and Teller, 2003; Rosentau et al., 2004; Teller and Leverington, 2004; Clarke et al., 2004; Rayburn et al., 2005; Leverington and Matile, 2006; Birks et al., 2007; Jakobsson et al., 2007; Makiaho, 2007; Clark et al., 2008).
Late Quaternary landscapes of isostatically-deformed regions can be reconstructed most easily and accurately through the employment of GIS-based techniques, since these techniques facilitate the collection of relevant topographic and isobase data, can be used to interpolate these values to continuous surfaces, and can be used to perform spatial arithmetic (Leverington et al., 2002b). In addition, quantitative analyses of resultant paleo-topographic databases (e.g., calculations of area, volume, aspect, and slope) can be performed using standard GIS tools, and specialized data products can be generated including topographic cross-sections, databases of hypsometry, determinations of drainage-network geometries, and three-dimensional visualizations. The raster data format is the most appropriate for use in the reconstruction of late Quaternary landscapes because it is ideally suited for arithmetic involving multiple data layers.
The use of GIS techniques in the reconstruction of a glacio-isostatically deformed landscape first involves the generation of a data layer for each of 1) modern topography; and 2) interpolated isobase values (or, alternatively, interpolated shoreline-elevation point data). If these data layers are georeferenced to each other, have been resampled to a common grid resolution, and contain values expressed with respect to a common datum, then a database of paleo-topography, in which elevations and bathymetry are given with respect to the restored water plane, can be generated by subtracting isobase values from modern elevations. This calculation adjusts topography for the effects of isostatic rebound since the time period being modelled, and also adjusts for other effects such as changes in sea level (Figures 1 and 2). Paleo-topographic databases, generated through the adjustment of modern elevations using strandline data, can be effectively used in the quantitative characterization of terrain morphology and surface processes (e.g., Lambeck, 1996; Gareau et al., 1998; Mann et al., 1999; Leverington et al., 2000, 2002ab; Schaetzl et al., 2002; Leverington and Teller, 2003; Rosentau et al., 2004; Teller and Leverington, 2004; Clarke et al., 2004; Yang and Teller, 2005; Rayburn et al., 2005; Leverington and Matile, 2006; Birks et al., 2007; Jakobsson et al., 2007; Makiaho, 2007; Clark et al., 2008; McCoy and Breckenridge, 2008).
Reconstructions based on field-derived isobase data are tied to real-world strandline features that define past water planes, providing a strong basis for inference of paleo-topography as long as strandline segments can be correlated and dated with sufficient confidence. Reconstruction results can alternatively be generated on the basis of geophysical principles of crustal loading and unloading (e.g., Tarasov and Peltier, 2006) if sufficient information is available regarding crustal properties and the nature of changes in ice loading through time, allowing for an independent and complementary means for arriving at analogous research outcomes. Importantly, for a given area of interest, recognition of the differences between the results of field-based and geophysics-based reconstructions has the potential to offer additional insights into the mechanisms of crustal rebound and associated environmental change.
Figure 1: Schematic representations of modern and ancient topography in a region influenced by late-Quaternary glacial loading (after Leverington et al., 2002b). A: Cross-section of modern topography and a deformed former water plane. B: Cross-section of paleo-topography of the same region after the effects of later crustal rebound have been removed.
Figure 2: Oblique views of schematic visualizations of the GIS method for reconstruction of paleo-topography from modern topography and databases of interpolated isobase data (after Leverington et al., 2002b).
Low-resolution global databases include the GLOBE database [modern land data; 30 x 30 arc-second cells (GLOBE Task Team, 1999)], the GTOPO 30 database [modern land data; 30 x 30 arc-second cells (EDC DAAC, 1996)] and the ETOPO5 and ETOPO2v2 databases [modern land and ocean data; 5 x 5 and 2 x 2 arc-minute cells (NGDC, 1988, 2008)]. These databases are useful in the reconstruction of the paleo-topography of continent-scale areas of interest, with pixel sizes generally no smaller than ~1 km. Databases with much higher spatial resolutions are also available, and allow for generation of paleo-topographic databases with superior representation of the local-scale variations of topography that are so critical to the proper identification and understanding of ancient shorelines and outlets. High-resolution databases of modern topography include those generated by NASA’s Shuttle Radar Topography Mission (SRTM) (e.g., 3 x 3 arc-second databases, with pixel sizes of ~90 m, available at no cost from institutions including the USGS; e.g., JPL, 1998), and topographic databases generated from contour-based topographic maps originally produced through photogrammetric analysis of aerial photographs (e.g., 3 x 3 arc-second databases available at no cost for the Canadian landmass from the Centre for Topographic Information, Natural Resources Canada; CTI, 2008).
Once entered into an appropriate GIS-readable data format from existing maps, isobase data (or, alternatively, shoreline-elevation point data) must be interpolated to a continuous surface and then resampled to the grid resolution of the database of modern elevations. The interpolation of isobase lines produces a continuous isobase surface, which can be used to describe the relative glacio-isostatic rebound that has occurred over a region since a given point in time. Within the geographical limits of the water body whose shorelines were used to generate the isobase map, the geometry of an isobase surface will correspond to the geometry of the corresponding (and now isostatically-deformed) water plane (Leverington et al., 2002b).
The Triangulated Irregular Network (TIN) algorithm (Peuker et al., 1978) has been generally found to produce good results when generating an isobase surface from point or curve data (Mann et al., 1999; Leverington et al., 2000, 2002b; Teller and Leverington, 2004), and Gaussian trend surfaces have in some cases been found to successfully fit rebound data to appropriate surfaces (Fretwell et al., 2004). Advanced geostatistical interpolation techniques such as kriging (e.g., Burrough and McDonnell, 1998) offer important alternatives to the TIN algorithm, especially when interpolations must be produced from relatively sparse isobase datasets. Geostatistical methods of interpolation can be effective in producing smooth isobase surfaces from widely-spaced data, whereas the suitability of the TIN model in such cases can be greatly impaired by the planar geometry of component triangular facets.
Subtraction of a raster database of interpolated isobase values from a database of modern topography will adjust modern topography for the effects of later glacio-isostatic deformation, allowing for the visualization and characterization of reconstructed landforms of interest.
A low-resolution demonstration of the generation of a database of paleo-topography using the GIS method described above is given here for the central Canadian Arctic (Figure 3) at 9300 14C yr BP (Leverington, 2002; Leverington et al., 2002b). At its maximum during the late Wisconsinan, the Laurentide Ice Sheet extended from the Arctic mainland of Canada across Victoria, Prince of Wales, and Somerset islands to the southernmost extents of Melville and Bathurst islands (e.g., Bryson et al., 1969; Dyke et al., 1982; Dyke, 1984; Hodgson et al., 1984; Dyke, 1987; Dyke and Prest, 1987; England et al., 2006). The Innuitian Ice Sheet extended across much of the Queen Elizabeth Islands during the late Wisconsinan (e.g., Mayewski et al., 1981; Dyke, 1984; England, 1998; Dyke, 1999; England, 1999; Wolfe and King, 1999; Dyke, 2003; England et al., 2006).
In order to reconstruct the low-resolution paleo-topography of the central Canadian Arctic at 9300 14C yr BP using the GIS method (Figure 3), an isobase surface and a database of modern topography were assembled for the region. Modern land elevations for the database of modern topography were extracted from the GLOBE database [cells in this database have a latitude-longitude grid spacing of 30 arc seconds (GLOBE Task Team, 1999)]. Modern bathymetric data were extracted from the ETOPO 5 database [cells in this database have a latitude-longitude grid spacing of 5 arc minutes (NGDC, 1988)]. The isobase surface was generated from the regional isobase map of Dyke et al. (1991) through digitization, TIN interpolation, and georeferencing of the interpolated database. The isobase surface, like the isobases from which it was derived, shows that the 9300 14C yr BP shoreline rises from low elevations of less than 50 m at Melville Island and eastern Devon Island, to elevations in excess of 250 m south of Boothia Peninsula.
Figure 3: A: Maps of the region used here for a low-resolution demonstration of the basic GIS-based technique for reconstruction of paleo-topography of deglacial landscapes. The approximate extent of glacial ice in the region at 9000 14C yr B.P. is given by blue lines (after Dyke and Prest, 1987; Dyke et al., 1996). B: A database of paleo-topography is generated through subtraction of an interpolated rebound surface (isobases depicted with white lines) from a database of modern topography (Leverington, 2002; Leverington et al., 2002b).
With both the database of modern topography and the isobase surface having identical geographic extents and cell sizes, a database of paleo-topography was produced by subtracting the isobase surface from the database of modern elevations. This database of paleo-topography delineates the approximate shorelines and subaerial extents of the region’s islands at 9300 14C yr BP, and is a topographic and bathymetric model of the region at that time. The depicted outline of sea level resembles that shown on the map for 9000 14C yr BP in Dyke and Prest (1987), and shows that much of the land area now exposed in the region was below sea level at 9300 14C yr BP. The surface areas of different land units in the region were affected in different ways by isostatic and sea-level influences at 9300 14C yr BP, depending on their topographic characteristics and the nature of local isostatic depression. While isostatic depression in the arctic islands greatly influenced the distribution of exposed land at 9300 14C yr BP, its associated influence on the morphology and size of major channels in the region was equally significant (Leverington, 2002; Leverington et al., 2002b).
As the southern margin of the Laurentide Ice Sheet (LIS) retreated during the last deglaciation of North America, meltwaters collected in proglacial lakes where drainage was impeded by ice. Lake Agassiz was likely the largest of these proglacial lakes, occupying at various times in its 4000-year history large portions of the southern part of the structural and topographic basin of central Canada. Because of its large area, volume, and overflow, Lake Agassiz had considerable influence on late-glacial North America (e.g., Teller, 1987; Kehew and Teller, 1994; Hu et al., 1997; Hostetler et al., 2000).
Reconstructions of the bathymetry of Lake Agassiz have been used to quantify aspects of its changing size and extent through time, and to estimate the volumes of outbursts driven by the deglaciation of new and lower outlets (Leverington et al., 2000, 2002; Teller and Leverington, 2004; Leverington and Matile, 2006) (e.g., Figures 6-11). A revision of aspects of the Lake Agassiz chronology depicted by Leverington et al. (2000, 2002a) and Teller and Leverington (2004) is due, in order to account for recent findings on the ground (e.g., Fisher, 2005, 2007; Fisher and Lovell, 2006; Rayburn and Teller, 2007).
Figure 4: Aerial photograph of preserved beach ridges near Miami, Manitoba, Canada.
Figure 5: Aerial photograph of Lake Agassiz spit formed in association with Bedford Hills Island, southeastern Manitoba, Canada.
Figure 6: Approximate total area occupied by Lake Agassiz over its 5000 calendar-year history (blue) (after Leverington and Teller, 2003). At no single time did the lake occupy this entire area; rather, the waters of the lake gradually shifted northward, following the wasting southern margin of the LIS and possibly combining with Lake Ojibway to the east (e.g., see Veillette, 1994; Leverington et al., 2002). Locations of lake overflow were a function of topography and the position of the confining ice sheet. Major drainage routes: S, Southern outlet to the Gulf of Mexico; K, Kaministikwia route to the Superior basin; NW, Northwestern outlet to the Arctic Ocean; E, Eastern outlet system to the Superior basin; KIN, Angliers and Kinojévis outlets to the St. Lawrence River valley; HB: final release of Lake Agassiz waters into the Tyrrell Sea.
Figure 7: Reconstruction of the bathymetry and shoreline of Lake Agassiz during the Herman stage (Leverington et al., 2000). This stage occurred at approximately ~11000 14C yr B.P. The volume of the lake at this time is estimated to have been about 10,900 km3, and the area of the lake is estimated to have been about 134,000 km2.
Figure 8: Reconstruction of the bathymetry and shoreline of Lake Agassiz during the Upper Campbell stage (Leverington et al., 2000). This stage occurred at approximately ~9400 14C yr B.P. The volume of the lake at this time is estimated to have been about 22,700 km3, and the area of the lake is estimated to have been about 260,000 km2. The dot near lake center marks the location of Winnipeg, Manitoba.
Figure 9: High-resolution reconstructions of paleo-topography allow for verification of flow routes inferred from field-based work. In the above example, segments of overflow routes connecting glacial Lake Agassiz (left) with glacial Lake Kelvin (right) are depicted (blue arrows) for the Upper Campbell stage of glacial Lake Agassiz (Leverington and Teller, 2003). Dashed lines represent positions of the southern margin of the Laurentide Ice Sheet in this region during the transition from the Upper Campbell stage to the Lower Campbell stage. Regions colored in green were above the elevation of Lake Agassiz water level during the Upper Campbell stage. Depicted flow routings are simplified versions of detailed drainage databases produced using the ArcHydro suite of tools.
Figure 10: High-resolution reconstruction of paleobathymetry of the Winnipeg region at the time of the Upper Campbell stage of glacial Lake Agassiz (Leverington and Matile, 2006). Terrain with elevations above the water level of Lake Agassiz at this time is colored in purple; the large purple zone in the southeastern part of the map is Bedford Hills Island. Overlain upon this reconstruction is a database of the field-mapped Lake Agassiz strandline features (red) of southeastern Manitoba, compiled from NATMAP geological databases.
Figure 11: Computer visualization of the Upper Campbell stage of Lake Agassiz, as viewed from 100 km above the future site of Winnipeg. The view is toward the northwest, with the southern extent of the Laurentide Ice Sheet forming the northern margin of the lake (compare this view with the map of the Upper Campbell stage given in Figure 8). The Manitoba Escarpment forms the green-shaded topographic high at bottom left. In the distance, the northwestern outlet of the lake is blocked by ice.
References Cited (c) David Leverington, 2009
(c) David Leverington, 2009