CC BY-NC-ND
0Contents lists available at ScienceDirect
Journal of Volcanology and Geothermal Research
journal homepage: www.elsevier.com/locate/jvolgeores
VOLCANOLOGY
Semi-automatic classification of glaciovolcanic landforms: An object-based mapping approach based on geomorphometry
G.B.M. Pedersen *
Nordic Volcanological Center, Institute of Earth Sciences, University of Iceland, Sturlugata 7,101 Reykjavik, Iceland
CrossMark
ARTICLE INFO
Article history:
Received 30 September 2015 Accepted 30 December 2015 Available online 15 January 2016
Keywords:
Glaciovolcanic landforms Tuya
Object-based image analysis (OBIA) Automated mapping techniques Geomorphometry Iceland
ABSTRACT
A new object-oriented approach is developed to classify glaciovolcanic landforms (Procedure A) and their landform elements boundaries (Procedure B). It utilizes the principle that glaciovolcanic edifices are geomorphometrically distinct from lava shields and plains (Pedersen and Grosse, 2014), and the approach is tested on data from Reykjanes Peninsula, Iceland. The outlined procedures utilize slope and profile curvature attribute maps (20 m/pixel) and the classified results are evaluated quantitatively through error matrix maps (Procedure A) and visual inspection (Procedure B). In procedure A, the highest obtained accuracy is 94.1%, but even simple mapping procedures provide good results (>90% accuracy). Successful classification of glaciovolcanic landform element boundaries (Procedure B) is also achieved and this technique has the potential to delineate the transition from intraglacial to subaerial volcanic activity in orthographic view. This object-oriented approach based on geomorphometry overcomes issues with vegetation cover, which has been typically problematic for classification schemes utilizing spectral data. Furthermore, it handles complex edifice outlines well and is easily incorporated into a GIS environment, where results can be edited or fused with other mapping results. The approach outlined here is designed to map glaciovolcanic edifices within the Icelandic neovolcanic zone but may also be applied to similar subaerial or submarine volcanic settings, where steep volcanic edifices are surrounded by flat plains.
© 2016 The Author. Published by Elsevier B.V. This is an open access article under the CC BY-NC-ND license
(http://creativecommons.org/licenses/by-nc-nd/4.0/).
1. Introduction
Maps of volcanoes represent a key information source for initial hazard assessment and they are therefore of major importance (e.g. Groppelli and Viereck-Goette, 2010). However, despite the wealth of remote sensing data available, mapping of volcanoes primarily relies on time-consuming, manual mapping. There is a considerable need for faster, automated mapping methods providing a low-cost route to volcanic mapping in remote or inaccessible environments. The objective of this study is to test two semi-automatic mapping methods applying object-based image analysis to slope and profile curvature layers derived from a medium-resolution digital elevation model (20 m/pixel). This is done in order to investigate the potential mapping of volcanic landforms based on their geomorphometric characteristics.
1.1. Volcanic mapping
The majority of volcanic mapping focuses on field-based mapping of the lithostratigraphy by observations of the lithology, the surface weathering characteristics, the outcrop patterns, and the degree of vegetation cover. These lithostratigraphic units can be defined more
» Tel.: +354 5254496; fax: +354 562 9767. E-mail address: grobirkefeldt@gmail.com.
comprehensively with various types of laboratory data, such as petro-graphic, geochemical, paleomagnetic, and radiometric analyses (e.g. Groppelli and Viereck-Goette, 2010). Furthermore, determination of the aerial extents and boundaries of the units are aided by aerial photointerpretation, which is particularly important where outcrops are sparse due to heavy vegetation or thick tephra cover (e.g. Neal and Lockwood, 2003; Herriott et al., 2008). Such geologic maps are of great importance and provide unprecedented detail of the volcanic deposits, allowing comprehensive characterization of eruptions and their timing. However, at the same time, this type of mapping is very costly and time consuming.
Satellite and airborne remote sensing (RS) have undergone a technical revolution providing a massive amount of data with high spatial, spectral, and temporal resolution via various sensors and satellite missions (e.g. Benediktsson et al., 2012). The wealth of RS data therefore provides an opportunity for detailed mapping, and Kervyn et al. (2007) demonstrated that diverse volcanic landscapes could be digitally mapped through visual interpretation of medium-resolution spectral and topographic satellite-based data. However, in order to exploit the full potential of the data, semi-automated- and automated mapping techniques are necessary.
1.2. Automated and semi-automated mapping of volcanic landforms
The land surface can be divided into a hierarchy of landscapes, land-forms, and landform elements, where the landform element is the
http://dx.doi.org/m1016/j .jvolgeores.2015.12.015
0377-0273/© 2016 The Author. Published by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).
smallest unit, indivisible at the given resolution and bounded by topographic discontinuities (Pike et al., 2009). Hence, volcanic landscapes can be divided into individual volcanic landforms (e.g. volcanic edifices), which in turn (dependent on the data resolution) can be subdivided further e.g. lava flows, lava channels, and levees.
Previous automated and semi-automated mapping techniques of volcanic landforms have focused on delineating the bases of volcanic edifices (Behn et al., 2004; Bohnenstiehl et al., 2012; Howell et al., 2012; Euillades et al., 2013). One approach has been the closed contour method, which locates quasi-circular topographic highs and subsequently selects the lowest, enclosing elevation (Behn et al., 2004; White et al., 2006; Bohnenstiehl et al., 2008; Cochran, 2008). The advantage of this method is that it operates directly on a digital elevation model (DEM), but a major drawback is that the base is constrained to have a constant elevation. This problem was overcome by Bohnenstiehl et al. (2012) who introduced adjustment of the base elevation along topographic profiles by the morphometric index cross-sectional area to outlining perimeter. However, the approach is still sensitive to user-defined parameters including the contour interval (Howell et al., 2012).
Other approaches use DEM-derived layers such as the slope, which is the first derivative, or the curvature. Curvature is a parameter that describes the concavity and convexity of a surface and is the second derivative of a DEM. The profile curvature is the curvature of the surface in the down-slope direction and affects the acceleration and deceleration of flow, and thereby influences erosion and deposition of material. Thus, positive profile curvature values indicate an upwardly concave surface, also called foot slopes. Negative profile curvature values denote shoulders, which are upwardly convex surfaces, and zero values indicate that the surface is linear (e.g. Fig. 8, Pedersen and Grosse, 2014a). On the other hand, planform curvature is the curvature of the surface perpendicular to the direction of maximum slope and influences convergence and divergence of flow. Hence, a negative planform curvature indicates a sidewardly concave surface, a positive planform curvature indicates a sidewardly convex surface, and zero denotes a linear surface.
The aim for methods using slope and curvature layers derived from DEMs is to delimit the base of volcanic edifices through changes in slope, since landform boundaries generally coincide with changes in slope and hence significant positive or negative curvature values (e.g. Dragut and Blaschke, 2006; Minar and Evans, 2008; Evans, 2012).
DEM-derived slope and profile curvature maps have been used for the identification of excavational volcanic landforms such as maars (Seib et al., 2013). Grosse et al. (2009, 2012) suggested a concave delimitation method for positive volcanic edifices by manual slope-break tracing on a combined slope-profile curvature map. This allows a systematic and uniform comparison of a variety of volcanic edifices in different geologic settings such as cinder cone fields in Mexico (Di Traglia et al., 2014) and glaciovolcanic edifices in Iceland (Pedersen and Grosse, 2014a). Euillades et al. (2013) automated this method and developed the NETVOLC algorithm, which automatically traces the concave edifice boundary by applying minimum cost flow networks. This method has been used for edifice delimitation in the near-global database on morphometry of composite volcanoes (Grosse et al., 2013) and has the advantage that it does not depend on several user-defined parameters. However, one of the caveats is that the algorithm requires input coordinates of the approximate center of the volcano.
Common to all these methods are that (1) they only focus on mapping the volcano edifice base and do not map individual volcano landform elements, (2) they only work for identification of isolated volcanic edifices, and (3) they yield erroneous results for complex edifices. This is a problem because many volcanic landscapes often consist of volcanic landforms with complex boundaries due to superimposed volcanic and tectonic structures.
The objective of this contribution is therefore to develop a mapping technique that addresses the abovementioned issues. This is achieved
by applying object-based image analysis (OBIA) to digital elevation model derived layers, such as slope and curvature maps. OBIA has the advantage that it can be applied to multiple scales and additional data can easily be added to the mapping procedure. Furthermore, the classified results can be imported directly to geographic information systems (GIS), and incorporated in a general mapping procedure along with other landforms.
2. Study area and data
The Reykjanes Peninsula, south-west Iceland (Fig. 1), is among the youngest and most pristine parts of Iceland and hosts a variety of well-preserved subaerial and glaciovolcanic edifices. The peninsula is primarily covered by basaltic lava flows that erupted after the termination of the last glaciation, estimated at around 12,000-15,000 years ago (e.g. Jakobsson et al., 1978; Ssmundsson et al., 2010). The glaciovolcanic edifices, on the other hand, were formed in contact with or confined by ice, resulting in distinct morphology and lithofacies (e.g. Noe-Nygaard, 1940; Matthews, 1947; Van Bemmelen and Rutten, 1955; Kjartansson, 1966; Jones, 1969). Most of these glaciovolcanic edifices are thought to be from either Early or Late Weichsel, although some deposits are older and have been ascribed to Early Brunhes (Ssmundsson et al., 2010). Previous geomorphometric analysis of basaltic volcanic edifices on the peninsula has shown that subaerial and glaciovolcanic edifices can be distinguished based on slope and profile curvature. This encouraged the investigation of whether a quantitative morphometric classification was possible (Fig. 8, Pedersen and Grosse, 2014a).
Recently, a geologic map of the peninsula was published by Ssmundsson et al. (2010) at a 1:100,000 scale and this map is used as a reference map to test the accuracy of the classification results. This geologic map was chosen because it covers the entire study area with the highest resolution. The map shows that >95% of the peninsula consists of two volcanic units: hyaloclastite and lava. The term hyaloclastite is used as a general term for hyalotuff, hyaloclastite, lapilli tuff, and pillow- and tuff-breccia. The lava is produced under subaerial eruption conditions, while the hyaloclastite is produced in subglacial, intraglacial, and submarine eruption conditions. The glaciovolcanic edifices can consist of only hyaloclastite (e.g. Sandfell, Fig. 2A-1), or both hyaloclastite and a lava cap (e.g. Geitafell, Fig. 2J-S). This depends on whether the eruption was purely intraglacial, or if the eruption protruded through the ice and produced a subaerial lava cap (for illustration, see Fig. 1, Pedersen and Grosse, 2014a). The glaciovolcanic edifices have diverse surface cover ranging from bedrock (either hyaloclastite or lava), to loose gravel and various types of vegetation (Fig. 2). This presents a significant problem when using spectral data for classification of these edifices (see Fig. 2 E-F and N-O).
The data used for this study are derived from a 20 m resolution DEM based on photogrammetry of aerial images spanning the time period from 1996 to 2012 (data were provided by Loftmyndir ehf). This resolution does not allow identification of individual lava flows or fissure swarms, but it is adequate for distinguishing the topographically distinct glaciovolcanic edifices down to ~0.1 km2 (i.e. 250 pixels). Hence, in this study, we distinguish between two classes: hyaloclastite (used as a general term) and lava fields (not distinguishing between individual flows).
3. Methods
3.1. Rationale
Classification of glaciovolcanic edifices using spectral data is, as mentioned, problematic due to the diverse edifice surface cover. However, Pedersen and Grosse (2014a) suggested that DEM-derived layers such as slope and profile curvature can be utilized. These authors defined a 5° gap in the average slope values (for a 20 m resolution
22°40'0nW 22°20,0"W 22WW 2r40'0"W
-1-1-1-1-
22°40'0"W 22°20'0"W 22*WW 2r40'0"W
Fig. 1. DEM-derived topographic map of the Reykjanes Peninsula ranging from 0 m (blue) to 500 m (red). The inset is a DEM of Iceland showing the location of the study area (red frame).
DEM) between the basaltic lava shields and basaltic glaciovolcanic edifices. Furthermore, they found that positive and negative profile curvature values within each glaciovolcanic edifice allow distinction of individual landform elements. Hence, negative profile curvature values coincide with the boundary between the flank and plateau, which for the lava-capped tuyas marks the boundary between intraglacial and subaerial volcanic activity in orthographic view (Fig. 8, Pedersen and Grosse, 2014a). Based on the abovementioned arguments, DEM-derived slope and profile curvature layers were therefore selected as input variables for the classification.
An object-based image analysis (OB1A) approach was chosen for classification rather than a pixel-based approach. OBIA is a two-step classification technique: (1) segmentation of an image creating segments/objects, each of which consists of a group of contiguous grid cells from a raster layer, followed by (2) rule-based classification of the individual objects (e.g. Benz et al., 2004).
This choice was made for a number of reasons; (A) OBIA classification allows incorporation of contextual information (e.g. texture, shape, directionality, spatial distribution within the study area, connectivity, etc.) and thereby opens up an opportunity to apply various geographical concepts to the image processing. Pixel-based classification, on the other hand, is only based on the pixel value and ignores spatial autocorrelation (Blaschke and Strobl, 2001; van Asselen and Seijmonsbergen, 2006; Dragu and Eisank, 2011; Van Den Eeckhaut et al., 2012). (B) A typical problem with pixel-based classification is that the classification results suffer from the salt-and-pepper effect due to scattered classification of single pixels. Since OBIA classifies objects rather than pixels, this problem becomes insignificant (e.g. Blaschke, 2010). (C) OBIA can, unlike pixel-based approaches, be applied on multiple scales and can thereby operate on objects of different sizes. (D) By using DEM-derived layers such as slope and profile curvature as input variables for the segmentation, the resulting objects are morphologically homogenous. The objects are thereby consistent with the landform elements defined by Speights (1974) or form facets as defined by Dikau (1989) (Romstad and Etzelmüller, 2012).
Automatic extraction of landforms and landform elements can be achieved by employing various terrain segmentation algorithms that
comprise either edge-based or region-based techniques. The edge-based techniques are preferable where significant surface discontinuities coincide with hydrological objects, while the region-based techniques aggregate pixels into groups according to user-prescribed rules of homogeneity. In this analysis, an eCognitions multi-resolution segmentation (Benz et al., 2004) has been used, which is a region-based technique that start with one-pixel objects.
Primarily, this type of OBIA-based classification has focused on fluvial, periglacial, and glacial landscapes (van Asselen and Seijmonsbergen, 2006; Romstad and Etzelmüller, 2009; Anders et al., 2011; Evans, 2012). It has also been used for landslide identification (Van Den Eeckhaut et al., 2012) as well as classification of topographic landforms at a global scale (Dragu^ and Eisank, 2012).
3.2. Procedure
The detection of landform and landform elements was conducted in eCognition 8.8 software (Trimble) (Benz et al., 2004). This provides an OBIA modular programming environment, where hierarchical rule sets of customized algorithms are built. Various rule sets were tested in order to explore the effect of parameter space and hierarchy, but overall, two main categories have been tested (Fig. 3):
(1) Slope-map-based procedure
(2) Combined profile curvature- and slope-map-based procedure
Procedure A extracts polygons with certain slope characteristics. It has the advantage that it only relies on the first derivative of the DEM, and therefore is less sensitive to DEM artifacts than procedure B, which incorporates the second derivative of the DEM (profile curvature). Ideally, it would be desirable to outline landforms and landform elements based on slope discontinuities rather than certain slope thresholds because it is discontinuity that defines the landform element boundaries. Therefore, a combined profile curvature and slope map procedure (Procedure B) was tested in order to extract polylines marking the boundaries between landform elements.
Fig. 2. Two examples of glaciovolcanic edifices: Sandfell (A-I) and Geitafell (J-S). (A&J) show a profile view of the two edifices. (B-F) and (K-Q) show maps of topography, slope, profile curvature, false-color multispectral image from SPOT satellite, and aerial image. The slope maps range from green (0°) to red (30°) and the profile curvature maps range from blue (—4) to red (4). The elevation scales for the topographic maps vary. Map view for Sandfell is 12 km across; map view for Geitafell is 36 km across. (G-I) and (P-S) show examples of the diversity of surface units in the field, and their location in the field can be seen on (B & K).
A slope map and a profile curvature map were created from the DEM in ArcGIS 10 using the Spatial Analyst module (Fig. 3). The slope map is the first derivative of the DEM and is calculated as the maximum rate of change in height for each cell with respect to a 3 * 3 cell neighborhood using the average maximum technique (Burrough and McDonnell, 1998). The profile curvature is the second derivative of a fourth-order
input surface fitted to a 3 * 3 window and was chosen since it affects the acceleration and deceleration of flow, and thereby influences erosion and deposition of material.
The initial step involves geometric terrain segmentation, where the slope (procedure A) or the profile curvature map and slope map (procedure B) serve as inputs for a multi-resolution segmentation,
Fig. 3. Generalized flowchart of procedure A and B. The blue boxes indicate that the module was performed in ArcGIS, while the yellow boxes indicate that the module was performed in eCognition. Used parameter space is mentioned in italic. O.L. is the abbreviation for object level.
MRS, creating the first object level (Fig. 4). The MRS in eCognition is a bottom-up region-growing algorithm where homogenous objects are formed starting with one-pixel objects that are subsequently merged into larger ones through a pair-wise clustering process (Benz et al., 2004). This clustering process stops when the object size exceeds the defined scale parameter and the clustering process minimizes the weighted heterogeneity of neighboring objects, taking slope and shape (Procedure A) or slope, profile curvature and shape (Procedure B) into account. The scale parameter is an abstract term and defines the maximum standard deviation of the homogeneity criteria and thereby determines the maximum allowed heterogeneity for the resulting image objects. Hence, at a given scale parameter, objects will be smaller for heterogeneous data than for more homogenous data.
Various scale and heterogeneity criteria (shape and compactness parameters) were tested through an error matrix (Procedure A, see Section 3.3) or adapted by visually inspecting the objects (Procedure
B). The objects were assigned classes based on slope thresholds, e.g. all objects with average slope > 12° are hyaloclastite (Procedure A) or profile curvature thresholds, e.g. all objects with average profile curvature >0.2 are concave (Procedure B) (Fig. 3). Small objects were subsequently removed and relative neighbor functions and pixel-based resizing (Module D-F in Fig. 3) were applied. Pixel-based resizing of objects allows generalization of object outlines by applying growing and shrinking algorithms that reshape the object outline by evaluating neighboring pixels according to a defined criterion (e.g. certain slope values). On the other hand, relative neighbor functions allow classification of neighboring objects which share a defined percentage of boundary with other classified objects. Furthermore, additional chessboard segmentation and subsequent reclassification with slope threshold values were implemented for procedure B, in order to obtain a better match with the reference data and improve the exported results (Fig. 3B, module G-H; Fig. 4, see paragraph 4.2 and Fig. 8). It is worth
Fig. 4. Overview of the hierarchal network showing the various levels in procedures A and B. The lower level (bottom) is the pixel level, which is the starting point for the MRS algorithm that generates object level 1. Subsequently, the objects are classified and merged (object level 2). In procedure B, this is followed by an additional chessboard segmentation, classification, and merged segmentation to improve classification results.
noting that modules D-F in procedure A and modules D-I in procedure B both use spatial information for classification (e.g. neighbor functions and object size) and multi-scale segmentation, which would not be feasible without an object-based approach.
Finally, the classified objects were exported either as polygon shapefiles (Procedure A) or as polyline shapefiles (Procedure B) (Benz et al., 2004) and imported directly to GIS for final editing.
3.3. Error matrix
To test the accuracy of procedure A, each classification result was tested against the reference map (Ssmundsson et al., 2010) by constructing an error matrix (e.g. Congalton, 1991; Stehman, 1997; Foody, 2002). The error matrix procedure was set up in ArcGIS 10 using the model builder module (Supp. Fig. 1).
Each classification result and the reference map were converted to a 20 m resolution grid file, where the pixels were either classified as hyaloclastite or lava. For the reference map, the hyaloclastite and lava field pixels where assigned values of 10 and 0, respectively. On the classified map, they were assigned to 1 and 0, respectively. The two maps were co-registered and summed, thereby constructing the error matrix map consisting of pixels with the values 0,1,10, or 11 (correct classified lava field = true negatives, incorrect classified lava field = false positives, incorrect classified hyaloclastite = false negatives and correct classified hyaloclastite = true positives). The proportion of correctly classified hyaloclastite and lava fields:
P(c) = P(0) + P(11)
gives the proportion of overall correctly classified areas and describes the overall map accuracy. The users accuracy, P(Ui,j), and producers accuracy, P(Pij), for classes i and j are defined as
P (Ui) = P(i)/P (i+)
P (Pi)= P( j)/P( j+)
k = (P(c)-YH1P(k+)* P(+k))/(1-^1=^+)* P(+k)),
where P(k+) = Yj^pkj and P(+k) =X==1pik
where P(i+) is the row marginal proportion and P(j + ) is the column marginal proportion. These accuracy measures have been chosen since they are fundamental and relevant for the direct interpretation ofthe data quality of each map (Stehman, 1997).
4. Results and discussion
The results are presented in two subsections, illustrating the results from procedures A and B, respectively.
4.1. Procedure A
To find the best rule set, various combinations of variables were evaluated, by changing one variable at a time (Fig. 3). For completeness, all variables in the MSR algorithm (scale, shape, compactness) were tested for a wide parameter range. Likewise, a fairly broad range of slope thresholds (used for classification) were tested as was the size for object removal, which heavily influences the number of incorrectly
classified small objects. In total, more than 700 different rule sets were evaluated individually (Supp. Table 1). The best classification result has an overall accuracy, P(c), of 94.08% and the corresponding error is displayed in Table 1, while the classification map is shown in Fig. 5.
However, very good overall accuracy (P(c) is> 90%) is obtained using the simple segmentation and classification procedure (without module E&F, Fig. 2), as shown in Fig. 6 (histogram of the overall accuracy for all the rule sets).
This is important because it signifies that fine-tuning of parameter spaces is not vital to obtain good results. The most influential parameters are slope threshold value and the object removal size. Generally, a slope threshold for object slope average > 12° and an object removal threshold of 500 pixels provide the best results, while a slope threshold for object slope average > 8° and an object removal threshold of 0 pixels display the poorest results (Fig. 7). The effects of compactness and shape parameters are negligible, whereas the scale parameter does have an effect, with optimal results for scale values of3,6, and 9 (Fig. 7).
Overall there is a good correspondence between the classified results and the reference map, showing that both fairly simple hyaloclastite edifices (Fig. 5 B-E) can be mapped, as well as edifices with very irregular and complex outlines (Fig. 5F-I). The user and producer accuracy of the lava plains is very high (95.85 and 97.74, respectively), while user and producer accuracy for hyaloclastite is moderately good to mediocre (70.47 and 56.04). The main reason for this is that this procedure relies on hyaloclastite units having steep slopes, which generally is a valid assumption (Pedersen and Grosse, 2014a). However, some hyaloclastite units that are older than Bruhnes have been eroded heavily by multiple glaciations, making them as flat as the surrounding lava plains (e.g. area around Engjahver, Mohalsdalur, Nupafjall, and Selfjall). It is also problematic when the hyaloclastites have been partly submerged by postglacial lavas, making the exposed hyaloclastite units topographically less distinct. Lava flows that have run down steep hyaloclastite slopes creating lava falls pose a different problem because they are steep and will be incorrectly classified as hyaloclastite (e.g. the lava falls in Fagridalur, Slatturdalur, and NatthagaskarS).
Finally, it is worth noting that the calculated error matrices assume that the reference map (Ssmundsson et al., 2010) is 100% correct, which is very likely not the case. A geological map at 1:100.000 scale is generalized and certain map units will have been merged in order to simplify the map. A likely consequence is that smaller-scale features and the exact location of unit boundaries may introduce an error in our error matrix.
4.2. Procedure B
Mapping landform element boundaries poses a greater challenge. This is because boundaries comprise a much smaller area than the land-forms themselves and because of the limitations of the export algorithm that convert the classified results into polylines. The classification procedure therefore has three additional modules (Fig. 3B) and involves 4 different object levels (Fig. 4B) in order to ensure good results.
Unlike procedure A, the direct calculation of an error matrix is not possible, since landform element boundaries are linear features, whereas the units on the reference map (Ssmundsson et al., 2010) are
Error matrix for the best classification result. The rule set has the following parameters: scale = 3; shape = 0.9; compactness = 0.4; slope threshold > 12°; object removal size < 500; classification of relative neighbors for objects with relative border = 1 and object size <500 pixels; looped pixel-based growth for slopes > 8°).
Classified data Reference data Lava plains Hyaloclastite Total User's accuracy [%] Producer's accuracy [%] Overall accuracy [%] KHAT
Lava plains 3,898,905 168,598 4,067,503 95.85 97.74 94.08 0.59
Hyaloclastite 90,077 214,919 304,996 70.47 56.04
Total 3,988,982 383,517 4,372,499
Fig. 5. The best classification result of procedure A. (A) Reykjanes Peninsula. Black frames mark the zoom boxes below, while the classification procedure and exact parameters are displayed in the top right corner. Frames (B-E) are 3.8 km across, while (F-I) are 7.6 km across. The yellow polygons are the classified hyaloclastite units, while the red line outlines the units from Ssmundsson et al. (2010). Rows 3 and 5 show the error matrix maps, where gray (value 0) marks correct classified lava fields, blue (value 1) denotes incorrect classified lava field, yellow (value 10) delimits incorrect classified hyaloclastite, and red (value 11) marks correct classified hyaloclastite.
31 1 27 30 17
90 00 90.25 90.50 90 75 91 00 91.25 91.50 91 75 92.00 92.25 92 50 92 75 93.00 93,25 93 50 93.75 94.00 94 25 94 50 94.75 95.00
Percentage (%)
Fig. 6. Histogram showing the proportion of overall accuracy, P(c), for all the rule sets tested.
polygons. The results have therefore been inspected visually by overlaying the boundaries on the reference map.
The best results were obtained with a scale factor = 1, the lowest value possible, because the boundaries are small. Like procedure A, the impact of shape and compactness in the MRS algorithm is minor on the classified objects, while unclassified objects with nearly no profile curvature are more influenced by the changing shape and compactness parameters (Sup. Fig. S2). Objects with average profile curvature values <— 0.2 were classified as convex and profile curvature values >0.2 were
classified as concave (Sup. Fig. S2). The boundaries were improved further using pixel-based resizing and employing relative neighbor relationships. Consequently, additional chessboard segmentation and slope thresholding were applied in order to confine the boundaries to slopes between 5° and 15°. This thinning of the landform element boundaries is necessary in order to improve the exported results because broad edifice boundaries often diverge around small obstacles, which causes problems when the classification results are exported as polylines. There are three shapefile export options in eCognition (polygon, main
! ! i i
Compactness
0,920 --
N Slope 5 8" ♦Slope £ 10° "Slope 2 12» ¿Slope S 14"
0.940 -
Object removal size
Fig. 7. The proportion of overall accuracy, P(c) as a function of the variables; scale (fixed parameters: shape 0.9; compact 0.4; object removal size <500); compactness (fixed parameters: scale 3; shape 0.9; object removal size <500); shape (fixed parameters: scale 3; compact 0.4; object removal size < 500), and object removal size (fixed parameters: scale 3; shape 0.9; compact 0.4). The tested slope threshold values were 8, 10, 12, 14.
Fig. 8. Example from the Fagradalsfjall area of the export shapefile options: (A) polygons, (B) main polyline, and (C) skeleton polyline. These results are superimposed on hyaloclastite units (yellow) from Ssmundsson et al. (2010). Red denotes concave boundaries and blue marks convex boundaries. (D) shows the result of the trimmed skeleton polylines.
polyline, skeleton polyline); their strengths and weakness are illustrated in Fig. 8.
Exported polygons delineate the boundaries nicely but are inconvenient since boundaries should preferably be marked as lines and not polygons. The main polylines do not allow branching, and classified boundaries are therefore excluded in the export process, which is a major problem (Fig. 8B, west part of Fagradalsfjall). The skeleton polylines match the classified areas well, but have a lot of small branches resulting from wide boundaries. However, most of these were removed by deleting small branches <300 m (using the trim line tool, ArcGIS editing toolbox) and satisfactory results are obtained (Fig. 8D).
The fairly simple mapping procedure is proficient in delineating concave and convex objects, which typically coincide with concave edifice bases and convex lines such as ridge crests or plateau edges (Fig. 9). For single mountains, the edifice bases are generally nicely delineated (Fig. 9 B-E), and even for complex edifices with very irregular outline and cross-cutting relationships (Fig. 9F-I) the method delivers fairly satisfactory results.
Thorbjorn (Fig. 9B), Geitahlid (Fig. 9D), Geitafell (Fig. 9E), Fagradalsfjall (Fig. 9F), Blafjall (Fig. 9H) as well as Storimeitill and Litlimeitill (Fig. 9I) show good correspondence between convex objects and the edges of plateaus. For lava-capped flat-topped glaciovolcanic edifices (e.g. Geitahlid, Geitafell, Fagradalsfjall, Blafjall,
and Storimeitill), these convex objects coincide with the mapped litho-logical boundaries between lava caps and hyaloclastites in orthographic view (S^mundsson et al., 2010), denoting the boundaries where the intraglacial volcanic activity changed to subaerial volcanic activity (Pedersen and Grosse, 2014a). In the case of Geitafell, multiple plateau levels are resolved, which had been simplified on the 1:100,000 map from S^mundsson et al. (2010) (see Fig. 2 J-S).
One problem with procedure B is that it classifies large faults and old sea cliffs. Furthermore, as for procedure A, problems arise when the hyaloclastite has very low slope values due to post-emplacement erosion (e.g. Keilir, Fig. 8C) or when younger lava has covered hyaloclastite slopes (e.g. Syllingafell in Fig. 8B).
5. Conclusions and perspectives
This study presents a new object-oriented approach to classifying glaciovolcanic landforms (Procedure A) and their landform elements boundaries (Procedure B). It is based on the principle that glaciovolcanic edifices are geomorphometrically distinct from the surrounding lava plains and is tested on the Reykjanes Peninsula, Iceland. The two procedures utilize slope (Procedure A) and a combination of slope and profile curvature attributes (Procedure B), which are derived from a 20 m/pixel resolution DEM.
Fig. 9. Classification results of concave (red) and convex (blue) boundaries. (A) Reykjanes peninsula. Black frames mark the zoom boxes below, while the classification procedure and exact parameters are displayed in the top right corner. Frames (B-E) are 3.8 km across; (F-I) are 7.6 km across. The yellow polygons are the mapped hyaloclastite units fromSsmundssonetal. (2010), which is displayed in rows 3 and 5. Postglacial lavas are marked with pink and purple colors.
Unlike previous automatic and semi-automatic methods employed to map volcanic landforms, this approach provides satisfactory results for very complex volcanic landforms and resolves sub-edifice scale structures. Procedure A was evaluated quantitatively by error matrixes and achieved accuracies as high as 94.1%, but even simple mapping procedures provide good results >90%. This demonstrates that fine-tuning of the method is not vital to obtain good results. The outputs from procedure B were evaluated by visual inspection and show positive results for classification of glaciovolcanic landform elements boundaries, including the boundary between hyaloclastite and lava caps.
The main limitation of both procedure A and B is that they do not classify topographically indistinct hyaloclastite, such as heavily eroded or partly submerged hyaloclastite. Despite this drawback, the major advantage of both methods is that they provide a reproducible, directly transferable method to map hyaloclastite units, glaciovolcanic edifice bases, and the edges of lava-capped plateaus on DEMs of similar resolution. The results are effortlessly implemented in a GIS mapping environment, and can therefore be easily edited and incorporated in a general mapping procedure along with other landforms.
Procedure A could be implemented in a broader mapping effort, to complete geologic mapping of the entirety of Iceland at 1:100,000 scale (only 30% of Iceland is presently covered (S^mundsson et al., 2010, 2012; Sigurgeirsson et al., 2015). Procedure B, on the other hand, provides an opportunity to uniformly map the lithological boundaries between lava caps and hyaloclastites in orthographic view throughout the Icelandic neovolcanic zone. The polylines can be converted to 3D polylines enabling transition heights from intraglacial to subaerial activity to be derived. The ability to map the spatial distribution of glaciovolcanic edifice volumes in multiple volcanic zones, as well as their transition heights, will also allow estimation of maximum jokullhlaup sizes for each of the analyzed glaciovolcanoes (Pedersen and Grosse, 2014b).
The mapping techniques outlined here are also of interest to the wider scientific community, for both planetary and submarine volcano research. Bathymetric DEMs are a primary data source for mapping and investigating submarine volcanoes, for determining both their distribution and erupted volumes. The glaciovolcanic edifices studied here show morphometric resemblances with submarine volcanoes (Gudmundsson and Jarosch, 2012), and hence the mapping techniques tested in this study are likely applicable to submarine volcanic settings as well.
Glaciovolcanic edifices have been proposed to exist on Mars (e.g. Chapman, 1994; Ghatan and Head, 2002; Pedersen et al., 2010; Pedersen, 2013), and given DEMs of comparable resolution, this method will provide an opportunity to make a uniform and direct comparison of the morphometric properties of terrestrial and potential Martian glaciovolcanic edifices.
Supplementary data to this article can be found online at http://dx. doi.org/10.1016/j.jvolgeores.2015.12.015.
Acknowledgments
The author would like to thank the Iceland Geosurvey, ISOR, for kindly providing the digital shapefiles from the geologic map of Ssmundsson et al. (2010), which has been essential for testing the classification results presented here. Furthermore, the author is very grateful for the constructive help and critique on the rule set development by Geoffrey Brian Groom from Biosciences, Aarhus University, Denmark, as well as very useful comments from reviewers, Lionel Wilson and Michelle Parks. Extended thanks go to Peder K. Bocher, Jens Christian Svenning, and Ecoinformatics at Biosciences, Aarhus University, Denmark, for input and inspiration on GIS modeling. This research was supported through grants from the Danish Research Council and the Nordic Volcanological Center, Institute of Earth Sciences, University of Iceland.
References
Anders, N.S., Seijmonsbergen, A.C., Bouten, W., 2011. Segmentation optimization and stratified object-based analysis for semi-automated geomorphological mapping. Remote Sens. Environ. 115 (12), 2976-2985.
Behn, M.D., Sinton, J.M., Detrick, R.S., 2004. Effect of the Galápagos hotspot on seafloor volcanism along the Galápagos Spreading Center (90.9-97.6°W). Earth Planet. Sci. Lett. 217 (3-4), 331-347.
Benediktsson, J.A., Chanussot, J., Moon, W.M., 2012. Very High-Resolution Remote Sensing: Challenges and Opportunities [Point of View]. Proceedings of the IEEE. 100 (6), 1907-1910.
Benz, U.C., Hofmann, P., Willhauck G., Lingenfelder, I., Heynen, M., 2004. Multi-resolution, object-oriented fuzzy analysis of remote sensing data for GIS-ready information. ISPRS J. Photogramm. Remote Sens. 58 (3-4), 239-258.
Blaschke, T., 2010. Object based image analysis for remote sensing. ISPRS J. Photogramm. Remote Sens. 65, 2-16.
Blaschke, T., Strobl, J., 2001. What's wrong with pixels? Some recent developments interfacing remote sensing and GIS. GIS-Z. Geoinformationssysteme 6, 12-17.
Bohnenstiehl, D.R., Howell, J.K., Hey, R.N., 2008. Distribution of axial lava domes along a superfast overlapping spreading center, 27-32°S on the East Pacific Rise. Geochem. Geophys. Geosyst. 9 (12), Q12016.
Bohnenstiehl, D.R., Howell, J.K., White, S.M., Hey, R.N., 2012. A modified basal outlining algorithm for identifying topographic highs from gridded elevation data, part 1: motivation and methods. Comput. Geosci. 49 (0), 308-314.
Burrough, P.A., McDonnell, RA, 1998. Principles of Geographic Information Systems. Oxford University Press, New York
Chapman, M.G., 1994. Evidence, age, and thickness of a frozen paleolake in Utopia Planitia. Mars. Icarus 109, 393-406.
Cochran, J.R., 2008. Seamount volcanism along the Gakkel Ridge, Arctic Ocean. Geophys. J. Int. 174(3), 1153-1173.
Congalton, R.G., 1991. A review of assessing the accuracy ofclassifications of remotely sensed data. Remote Sens. Environ. 37,35-46.
Di Traglia, F., Morelli, S., Casagli, N., Monroy, V.H.G., 2014. Semi-automatic delimitation of volcanic edifice boundaries: validation and application to cinder cones of the Tancitaro-Nueve Italia region (Michoacán-Guanajuato Volcanic Field, Mexico). Geomorphology 219,152-160.
Dikau, R., 1989. The application of a digital relief model to landform analysis in geomorphology. In: Raper, J. (Ed.), Three Dimensional Applications in Geographical Information Systems. Taylor and Francis, London, pp. 51 -77.
Draguf, L., Blaschke, T., 2006. Automated classification of landform elements using object-based image analysis. Geomorphology 81 (3-4), 330-344.
Draguf, L., Eisank, C., 2011. Object representations at multiple scales from digital elevation models. Geomorphology 129 (3-4), 183-189.
Draguf, L., Eisanlt, C., 2012. Automated object-based classification of topography from SRTM data. Geomorphology 141-142 (0), 21-33.
Euillades, L.D., Grosse, P., Euillades, P.A., 2013. NETVOLC: an algorithm for automatic delimitation of volcano edifice boundaries using DEMs. Comput. Geosci. 56 (0), 151-160.
Evans, I.S., 2012. Geomorphometry and landform mapping: what is a landform? Geomorphology 137 (1), 94-106.
Foody, G.M., 2002. Status of land cover classification accuracy assessment. Remote Sens. Environ. 180,185-201.
Ghatan, G.J., Head III, J.W., 2002. Candidate subglacial volcanoes in the south polar region of Mars: morphology, morphometry, and eruption conditions. J. Geophys. Res. 107 (E7), 5048.
Groppelli, G., Viereck-Goette, L., 2010. Introduction. In: Groppelli, G., Viereck-Goette, L. (Eds.), Stratigraphy and Geology of Volcanic Areas: Geological Society of America Special Paper http://dx.doi.org/10.1130/2010.2464(00).
Grosse, P., van Wyk de Vries, B., Petrinovic, I.A., Euillades, P.A., Alvarado, G.E., 2009. Morphometry and evolution of arc volcanoes. Geology 37 (7), 651-654.
Grosse, P., van Wyk de Vries, B., Euillades, P.A., Kervyn, M., Petrinovic, IA, 2012. Systematic morphometric characterization of volcanic edifices using digital elevation models. Geomorphology 136 (1), 114-131.
Grosse, P., Euillades, P., Euillades, L., van Wyk, D.V., 2013. A global database of composite volcano morphometry. Bull. Volcanol. 76 (1), 1-16.
Gudmundsson, M.T., Jarosch, A.H., 2012. Subglacial and subaqueous monogenetic volcanoes, thermal and morphological contrasts. AGU AnnualFall Meeting, V21A-2761.
Herriott, T.M., Sherrod, D.R., Pallister, J.S., Vallance, J.W., 2008. Photogeologic maps of the 2004-2005 Mount St. Helens eruption. In: Sherrod, D.R., Scott, W.E., Stauffer, P.H. (Eds.), A Volcano Rekindled: The Renewed Eruption of Mount St. Helens, 2004-2006. U.S. Geological Survey, pp. 209-224.
Howell, J.K., White, S.M., Bohnenstiehl, D.R., 2012. A modified basal outlining algorithm for identifying topographic highs in gridded elevation data, part 2: application to Springerville Volcanic Field. Comput. Geosci. 49 (0), 315-322.
Jakobsson, S.P., Jónsson, J., Shido, F., 1978. Petrology of the western Reykjanes peninsula, Iceland. J. Petrol. 19 (4), 669-705.
Jones, J.G., 1969. Intraglacial volcanoes of the Laugarvatn region, south west Iceland 1. Q. J.Geol. Soc. 124,197-211.
Kervyn, M., Kervyn, F., Goossens, R., Rowland, S.K., Ernst, G.G.J., 2007. Mapping volcanic terrain using high-resolution and 3D satellite remote sensing. Geol. Soc. Lond., Spec. Publ. 283 (1), 5-30.
Kjartansson, G., 1966. A comparison of table mountains in Iceland and the volcanic Island of Surtsey off the south coast of Iceland (english summary). Natturufroeingurinn 36, 1-34.
Matthews, W.H., 1947. "Tuyas", flat-topped volcanoes in northern British Colombia. Am. J. Sci. 245, 560-570.
Minär, J., Evans, I.S., 2008. Elementary forms for land surface segmentation: The theoretical basis of terrain analysis and geomorphological mapping. Geomorphology. 95 (34), 236-259.
Neal, CA, Lockwood, J.P., 2003. Geologic map of the summit region of Kilauea Volcano, Hawaii. USGS map Ml-2759.
Noe-Nygaard, A., 1940. Sub-glacial volcanic activity in ancient and recent times (Studies in the palagonite-system of Iceland no. 1) Tom. 1, no. 2 Folia Geogr. Dan. 1,1-67.
Pedersen, G.B.M., 2013. Frozen Martian lahars? Evaluation of morphology, degradation and geologic development in the Utopia-Elysium transition zone. Planet. Space Sci. 85 (0), 59-77.
Pedersen, G.B.M., Grosse, P., 2014a. Morphometry of subaerial shield volcanoes and glaciovolcanoes from Reykjanes Peninsula, Iceland: effects of eruption environment. J. Volcanol. Geotherm. Res. 282 (0), 115-133.
Pedersen, G.B.M., Grosse, P., 2014b. Insights into subglacial eruptions based on geomorphometry: broad scale analysis of subglacial edifices in Iceland. EGU General Assembly. p. 16.
Pedersen, G.B.M., Head Ill, J.W., Wilson, L., 2010. Formation, erosion and exposure of early Amazonian dikes, dike swarms and possible subglacial eruptions in the Elysium Rise/ Utopia Basin Region, Mars. Earth Planet. Sci. Lett. 294,424-439. http://dx.doi.org/10. 1016/j.epsl.2009.08.010.
Pike, R.J., Evans, I.S., Hengl, T., 2009. Geomorphometry: a brief guide. In: Hengl, T., Reuter, H.l. (Eds.), Geomorphometry: Concepts, Software, Applications. Developments in Soil Science vol. 33. Elsevier, pp. 3-30.
Romstad, B., Etzelmüller, B., 2009. Structuring the digital elevation model into landform elements through watershed segementation of curvature. Proceedings of Geomorphometry, p. 55.
Romstad, B., Etzelmüller, B., 2012. Mean-curvature watersheds: a simple method for segmentation of a digital elevation model into terrain units. Geomorphology 139-140 (0), 293-302.
Ssmundsson, K., Johannesson, H., Hjartarson, Ä., Kristinsson, S.G., Sigurgeirsson, MA, 2010. Geological map of Southewest Iceland. Iceland Geosurvey. 1:100000.
Ssmundsson, K., Hjartarson, Ä., Kaldal, I., Sigurgeirsson, M.Ä., Kristinsson, S.G., Vikingsson, S., 2012. Geological Map of the Northern Volcanic zone, Iceland. Northern Part. 1:100.000 (Reykjavik), Iceland GeoSurvey and Landsvirkjun.
Seib, N., Kley, J., Büchel, G., 2013. Identification of maars and similar volcanic landforms in the West Eifel Volcanic Field through image processing of DTM data: efficiency of different methods depending on preservation state. Int. J. Earth Sci. 102 (3), 875-901.
Sigurgeirsson, M. Ä., Hjartarson, Ä., Kaldal, I., Ssmundsson, K., Kristinsson, S. G., Vikingsson, S., 2015. Geological map of the Northern Volcanic Zone, Iceland. Southern Part. 1:100.000 (Reykjavik), lceland GeoSurvey and Landsvirkjun.
Speights, J.G., 1974. A parametric approach to landform regions. ln: Brown, E.H., Waters, R.S. (Eds.), Progress in Geomorphology: Papers in Honour of David L. Linton. Institute of British Geographers. Alden & Mowbray Ltd at the Alden Press, Oxford, pp. 213-230.
Stehman, S.V., 1997. Selecting and interpreting measures of thematic classification accuracy. Remote Sens. Environ. 62 (1), 77-89.
van Asselen, S., Seijmonsbergen, A.C., 2006. Expert-driven semi-automated geomorpho-logical mapping for a mountainous area using a laser DTM. Geomorphology 78 (3-4), 309-320.
Van Bemmelen, R.W., Rutten, M.G., 1955. Table mountains of Northern Iceland. E.J. Brill, Leiden, Netherlands.
Van Den Eeckhaut, M., Kerle, N., Poesen, J., Herväs, J., 2012. Object-oriented identification of forested landslides with derivatives of single pulse LiDAR data. Geomorphology 173-174 (0), 30-42.
White, S.M., Haymon, R.M., Carbotte, S., 2006. A new view of ridge segmentation and near-axis volcanism at the East Pacific Rise, 8°-12°N, from EM300 multibeam bathymetry. Geochem. Geophys. Geosyst. 7 (12), Q12005.
