Detection of crevassed areas with minimum geometric information: Vernagtferner case study

Abstract

Crevasses pose severe risks for mountaineers and field glaciologists. Smaller cracks between 0.5 and 2 m are still dangerous, but often not visible in medium resolution satellite imagery. If they are snow covered, they are completely undetectable by optical sensors. We set out to develop an approach to detect potentially crevassed areas by a minimum of geometric data, and to make the method generally applicable to glacier regions. On Vernagtferner, we compared a reference dataset of crevasses observed in high-resolution optical imagery with the curvature of the ice surface and the spatial gradients in driving stress. Both parameters can be derived from a digital surface model and a bedrock model, derived from ice thickness measurements. The correlation patterns show that crevasses preferably form in convex areas and in areas where the driving stress rapidly increases. This corresponds with the theory of crevasse formation. Although the method still misclassifies larger parts, the approach has the potential to define probable non-crevassed areas as well as to aid the planning of safe routes.

1. Introduction

Research on crevasses has a long tradition. Colgan and others (2016) provides an extensive overview of the most important previous studies. There are several reasons for the research efforts on crevasse distribution. Firstly, crevasses are hazardous for mountaineers. This makes it useful to know the location of crevasses with a minimum width of around 30 cm, large enough to cause a serious accident. Secondly, crevasse distribution influences various processes on a glacier. Crevasses can affect the glacier mass balance in different ways, e.g., crevasses lead to the increased collection of solar energy (Cathles and others, 2011). In addition, the presence of crevasses can be an indicator of surge dynamic process, which is an important, but little understood process in glacial acceleration (Trantow and Herzfeld, 2018). Englacial hydrology is also affected by crevasses and thus water transport to the glacier bed. This, possibly, has consequences for basal sliding and might be crucial to the surge behavior of certain glaciers (Fountain and Walder, 1998).

Several well-known methods are used to detect visible crevasses by airborne or spaceborne imagery, or indeed laser altimetry. The detectability of crevasses here depends on the resolution of the imagery or the spatial scanning rate of laser altimetry data. With the use of standard high resolution satellite data (1–10 m resolution), only relatively broad crevasses can be detected. The detection of crevasses from satellite data can be improved by using deep neural network (Lai and others, 2020). Currently, image resolutions are getting higher, making smaller crevasses detectable. In addition, higher resolutions can be achieved by airborne or unpiloted aerial vehicle (UAV) -borne images. However, collecting high-resolution data requires a significant amount of effort, e.g., by using drone flights for a large area. Buried crevasses can be located using active sensors in the micro wave band, where the wavelength limits penetration depth and spatial resolution. Nevertheless, this does not limit crevasse detection, as Marsh and others (2021) shows, it is possible to detect crevasses narrower than the spatial resolution or at a few meters depth, depending on the imaging parameters, such as water content of the snowpack or view direction of the sensor. In order to find small buried crevasses, ground penetrating radar (GPR) can also be used, which is very labor intensive and additionally carries the risk of a fall, in the case of surface-bound measurements. Therefore, this method is not feasible for large areas. Thompson and others (2020) analyzed satellite- and helicopter-based methods to separate a polar glacier into crevassed and non-crevassed areas. They come to the result, that the detection works with some limitations on the speed, the area and the crevasse geometry, such as crevasse width.

Apart from the attempts to detect existing visible or buried crevasses with sensors in the field, certain well-established physical approaches to modeling crevasses exist. According to Colgan and others (2016), these include, for example, the classical zero stress model as one of the oldest models, the fracture mechanics model and the continuum damage mechanics model. However, these models are affected by uncertainties and they require a considerable amount of specific input data, which must be very precise. Such information is usually not available for larger areas.

The zero stress model can be used to calculate the crevasse spacing from the strain rate, which can be derived by velocity data, so that a critical yield strain rate is employed (Clason and others, 2012). The fracture mechanics model is a parameter intensive approach, involving a significant amount of a priori information, e.g., the ice fracture toughness must be well-known (Clason and others, 2012). Iceberg calving and failure of hanging glaciers can be modeled by the continuum damage mechanics models, based on parameters determined in the laboratory (Pralong, 2005).

Nevertheless, we aimed to find a method to detect crevassed glacier areas that requires less effort and input data, hence also making the approach feasible for larger regions. In our research, we were successful in finding a parameter dependence for the location of crevassed regions using only a digital surface model (DSM) and in addition, an estimate of the ice thickness. From these data, the parameters curvature and change in driving stress (Δτ _di) were derived and their dependence on crevassed regions analyzed.

Previous research has already shown that crevasse formation does not depend on the existence of old crevasses (Whillans and others, 1993) and ‘ogives’ (Goodsell and others, 2002). In our approach to investigate crevasse formation, we use the geometry of the glacier, in contrast to most other studies that use the dynamics of a glacier and thus mostly use the critical strain rate or principal stress to infer crevasse formation (Ambach, 1968; Harper and others, 1998; van der Veen, 1999; Kaluzienski and others, 2019). Vaughan (1993) compared crevasse formation with measured strain rates by converting the strain rates to stresses. Using Glen's flow law, Vaughan showed that crevassing is dependent on the creep, but he concluded that the determination of a failure criterion is difficult. The establishment of a failure criterion has also been part of many other studies. For example, Koziol and others (2017) employed an ice yield strength of 132.5 kPa to determine crevassed areas after using velocities to derive surface stresses. In addition, Mottram and Benn (2009) used a critical strain rate to determine the crevasse depth, which was calculated experimentally and corresponds to a stress of about 60 kPa. Crevasse interaction and the crevasse spacing/depth ratio may be influenced by e.g. stress gradients, surface curvature, critical stress and other parameters (Sassolas and others, 1996). Other methods, such as the use of high resolution digital elevation model (DEM) derived from UAV have been done. Florinsky and Bliakharskii (2019) derived morphometric variables such as horizontal and vertical curvature from a high resolution DEM. In this case, crevasses must be visible or covered by very small snowpack so that the contours of the crevasses can still be seen in the DEM. The method is only suitable for crevasses with a minimum width of 1.5 m. Ryan and others (2015) used a slightly different method to detect crevasses from a high resolution DEM. They subtracted a filtered DEM from the original DEM to obtain negative surface anomalies, i.e. crevasses. The detection of the crevasses is strongly dependent on the resolution and accuracy of the DEM.

The driving stress is linked to other important glacier parameters, e.g., the driving stress has a significant connection with the acceleration (Scott and others, 2009). Hedfors and others (2003) examined the ratio and thus the interaction of basal drag and driving stress in context with the bed topography. Thus far, driving stress has not been in the focus of crevasse research, even though it is closely related to strain rates and might, therefore, also relate to fracture conditions.

2. Study area and data

2.1 Study site

For our investigations we selected Vernagtferner in Austria. It is located in the southern Oetztal Alps and extends over an elevation of ca. 2800 to 3500 m a.s.l.. Between 1981 and 2018 a cumulative mass loss of more than 26 000 mm w.e. was observed, which is closely related to the changing climate conditions in the region (BAdW, 2021). Despite this fact, it is still one of the largest glaciers in Austria, with a size of approximately 6.9 km². Due to its size and its extent across several undulated basins, Vernagtferner shows several distinct crevassed areas.

Vernagtferner has a long history of geodetic and glaciological surveys. Indeed, Sebastian Finsterwalder (1897) generated the first detailed map of this glacier between 1888 and 1889. Since then, a large number of research activities have concentrated on Vernagtferner, specifically since 1962, when the Commission for Glaciology (today: Geodesy and Glaciology) of the Bavarian Academy of Sciences and Humanities in Munich was founded (Braun and others, 2013). Bi-annual glacier mass balances have been determined since 1964 and run-off from the drainage basins has been measured since 1973.

This makes Vernagtferner one of the few glaciers worldwide where the hydrological balance can be compared with the glaciological balance over several decades. Apart from in-situ measurements, this glacier has been surveyed frequently by terrestrial and aerial photogrammetry, as well as laser scanning. This means high-resolution imagery is available for different time periods (Mayer and others, 2013).

2.2 Dataset

For our study, we used a DSM from an aerial photogrammetry survey in 2018 with a raster resolution of 1 m. In addition, the corresponding orthoimages (resolution 20 cm) were used to detect visible crevasses manually. The photogrammetric survey took place during minimal snow cover thickness and extent at the end of the ablation season. In order to analyze the physical conditions of the glacier, we also required information on the bedrock elevation. The respective data set results from differencing an ice thickness model of Vernagtferner from the DSM in the according year (2006).

The ice thickness map was derived from observed GPR profiles with a 40  MHz antenna during February 2006 (Mayer and others, 2013). The horizontal resolution was about 5 cm, however, the horizontal accuracy was only 4 m. The vertical accuracy was around 7 m and the vertical resolution was about 2 m. The spacing between the survey grid lines was variable, the mean spacing was about 100 m or less. The interpolation between the GPR profiles was carried out by a Kriging algorithm, yielding an overall accuracy of ±7 m. Some regions below Hochvernagtspitze, Hinterer and Vorderer Brochkogel were too difficult to access and could not be measured. These data gaps were interpolated, based on surface slope and mass conservation considerations. In general, the thickness in the non-observed parts was less accurate. This was also the case for regions close to the margins, where ice thickness had been manually set to zero for the interpolation process. The corresponding bedrock map served as a basis for deriving the ice thickness distributions for 2018 based on the existing DSM.

The derived ice thickness distribution from the bedrock DEM 2006 and the DSM 2018 is depicted in Fig. 1. Glacier retreat between 2006 and 2018 results in negative ice thickness at the present glacier margin and regions not coved by thickness measurements, because of interpolation issues of the 2006 bedrock model. The GPR grid for the bedrock model was designed to cover as much of the major ice body and its deepest parts, but not the marginal areas. The GPR device could only reliably detect ice thicknesses of more than approximately 15 m due to the low frequency used. The voids along the margins resulted in an underestimation of the ice thickness due to a rather linear interpolation from the GPR profiles towards the glacier margin in 2006. However, areas with high uncertainties in ice thickness are mainly limited to the peripheral zones, where fewer, or no, crevasses are located. Thus, the impact on this study can be disregarded. The method was applied only in areas with valid ice thicknesses.

Fig. 1. Overview of Vernagtferner in 2018 with its individual glacier tongues and additional local designations. The background orthoimage is based on an aerial survey in 2018.

2.3 Reference dataset of crevasses

The visible crevasses with a minimum width of ca. 0.3 cm could be identified on orthoimages with a spatial resolution of 20 cm. They were delineated manually for the year 2018 (Fig. 2a).

Fig. 2. (a) Manually detected crevasses as polygons in 2018 with a background orthoimage of 2018. (b) Crevasse index and crevassed areas considered in this study. Some areas were excluded due to unreliable ice thickness information.

Difficulties in the delineation persist due to partly snow-covered crevasses, which are located mainly along the upper glacier margin. Moreover, varying lighting conditions owing to shadows affect the identification of crevasses in some parts of the glacier.

The mapped crevasses were saved as polygons in a shape-dataset. Approximately 9 500 crevasses were detected for the year 2018. The high spatial resolution of 20 cm of the orthoimage results in a rather high accuracy of the crevassed area, including crevasse geometries.

The aim of this study was to investigate whether regions with high crevasse densities can be detected with a minimum of geometric information. For this purpose, the percentage of crevassed surface area within a preassigned 20 m × 20 m raster was used as an index of crevasse density.

In order to obtain a generalized crevassed density for each raster field and the neighboring region, a linear distance filter was applied on the calculated index. A weight-matrix was used as linear distance filter (Fig. 3) to consider the percentage of crevassed surface of the neighboring pixels with a decreasing weight with distance. The calculation of the crevasse index CI was done as presented in equation (1), with W representing the weight-matrix and A containing the crevassed surface in percent for each 20 m × 20 m raster.

(1)$$CI_{( i, j) } = \left({{\sum\nolimits_{k = -2}^{2} \quad \sum\nolimits_{l = -2}^{2} ( W_{( i-k; j-l) } \ \, A_{( i-k; j-l) }}) \over \sum\nolimits_{k = -2}^{2} \sum\nolimits_{l = -2}^{2} W_{( i-k; j-l) } } \right)$$

Fig. 3. Weight-matrix W for the calculation of a crevasse index CI for a 20 m x 20 m raster by applying a linear distance filter on the percentage of crevassed surface of the raster and its surrounding area.

As a result, the calculated index for each 20 m × 20 m raster includes its percentage of crevassed surface, as well as the surrounding cells in a 100 m × 100 m area. Figure 2b shows the calculated crevasse index for Vernagtferner. The crevasse index ranges from 0 (not crevassed) to 17 (highly crevassed). Finally, a threshold for the crevasse index was established to determine whether the 20 m × 20 m raster is an area with high crevasse densities or an area with less or no crevasses. The crevasse index was required to be at least 4 or higher to be classified as a crevassed area. Thus, 4% of the area in the wider region (100 m × 100 m area) needed to consist of identified crevasses after applying the linear distance filter to be considered a crevassed area. This crevasse index map was used as reference dataset in the subsequent analysis, excluding areas with unreliable ice thickness information due to their proximity to the glacier margin.

3. Methodology

We decided to investigate the relation of two parameters to the crevasse index, the spatial gradient in driving stress and the curvature. Both can be derived from a DSM, which is available in a sufficient resolution for many parts of the Alps. For the analysis, a high-resolution DSM with a raster resolution of 1 m was used, since its measurement date was similar to the reference dataset. However, it was downscaled to 20 m. Thus, other global DEMs, such as the Copernicus DEM (Copernicus, 2022) with a raster resolution of 25 m or the regional DEM Austria 10 m resolution (Austria, 2022) could also be used to apply our approach. The use of data with a higher resolution is neither of advantage for the driving stress nor the curvature, because these parameters are not effective at a very small scale. However, for the driving stress calculation, information on ice thickness was required in addition. This information is often not available from direct measurements, but there exist an increasing number of ice thickness estimates from different model approaches. Farinotti and others (2017) provides an overview of the previous approaches for the ice thickness estimation, using the characteristics of the glacier surface. Ice thickness information from a database, such as GlaThiDa (Welty and others, 2020) can also be used if it contains data with appropriate date and spatial resolution.

3.1 Driving stress

In response to the driving stress, glacier ice is deformed and ice transport is enabled. In general, this transport acts in the direction of the surface slope. The driving stress (τ _di) is a function of the bulk glacier density ρ, the gravity g, the ice thickness h and the glacier surface slope α:

(2)$$\tau_{di} = \rho \;\ast\; g \;\ast\; h \;\ast\; \sin\! \alpha$$

The main source of ice deformation is the driving stress, i.e. the residual shear stress component along the slope based on the general stress tensor. This is especially true for cold-based glaciers, where basal sliding can be neglected. The driving stress acts as a shear stress throughout the ice column and is counteracted mainly by basal drag at the glacier bed (apart from lateral and bridging effects within the ice column).

In order to relate the stress conditions to ice deformation, the driving stress needs to be connected to strain rates by a constitutive relation. For glaciers, Glen's flow law (Glen J, 1955) is usually used, where the strain rate is non-linearly dependent on the stress.

Crevasses occur when a critical tensile stress is exceeded. Due to the correlation between stress and strain, a critical strain rate for the formation of crevasses can be expected. The strain rate, in turn, is influenced by the creep, which results, among other things, from the driving stress (van der Veen, 1998; Colgan and others, 2016).

Therefore, a strong change in driving stress across short spatial scales might be connected to the occurrence of crevasses and might be a suitable indicator for crevasse formation. Moreover, we assume that compressive conditions are usually not connected to the generation of crevasses. Such conditions exist when the driving stress decreases along flow. We assume that ice flows in the direction of the steepest slope. Consequently, ice flow direction can be derived from a DSM and compared to the stress gradient in the same direction. Compressive flow conditions, i.e. a positive driving stress gradient along flow, resulted in neglecting these regions from the analysis. Therefore, we only consider regions where the driving stress generally leads to dilation of the ice.

We calculated the driving stress (equation (2)) based on the surface slope and ice thickness fields derived from the input data, using an ice density of 917 kg m⁻³ and a value of 9.81 m s⁻² for the gravitational acceleration. Thereby, the uncertainty of the driving stress is about 17.5 Pa. The driving stress was calculated for an 80 m × 80 m grid because the parameter is not effective at a smaller scale. For this purpose, Median and Gaussian filters were used to smooth the DSM and ice thickness during the resampling and remove possible outliers. The spatial change in driving stress, which may be relevant for the crevasse formation, was calculated by the hypotenuse of the gradients in easting and northing of the determined driving stress field.

3.2 Curvature

Curvature of the glacier surface is the second parameter whose influence on the formation of the crevasse fields is investigated. The curvature can be calculated by the second spatial derivative of the surface, while the first derivative represents the slope.

The parameter curvature has been used in previous studies to calculate the flow type of a glacier, which is relevant for crevassing. While a concave glacier surface results in compressive flow, a convex glacier surface supports extending flow (Nye, 1952). We assume that a strong surface curvature is related to strong spatial changes in the general stress conditions. This again might be connected to a potential increase of strain rates across the threshold of brittle failure and crack formation in the ice column.

In topographic analyses, a vertical curvature can be distinguished from a horizontal curvature. The first describes the curvature along the steepest slope, while the second describes the curvature perpendicular to the steepest slope. It is possible to determine a combination of vertical and horizontal gradients, which can represent many geological processes in alpine geosystems, as shown in other research studies (Rasemann, 2003; Dikau and others, 2019).

As the flow direction of a glacier also follows the steepest slope, we focus on the vertical curvature. Different definitions of vertical surface curvatures exist, which are compiled in Minár and others (2020). The vertical curvature tool of ArcGIS (version 10.8.1) was used for the analysis. Here, no geometric curvature was used, but the 2^nd slope line derivative (z _ss), which is defined in Minár and others (2020):

(3)$$z_{ss} = - {z_{xx}z_x^2 + 2z_{xy}z_xz_y + z_{yy}z_y^2\over ( z_x^2 + z_y^2) }$$

The 2^nd slope line derivative is a sub-form of the basic curvatures. Differences between sub-form and basic curvatures are eliminated due to the later categorization into the classes concave, flat and convex (Minár and others, 2020).

The vertical curvature was calculated from the 2018 DSM. Similar to the procedure used for the drive stress calculation, the DSM was resampled to 20 m × 20 m.

4. Results

4.1 Spatial change in driving stress

The spatial change of the driving stress might be indicative of crevasse formation if the spatial gradients in driving stress are strong. The gradients of the driving stress are presented in Fig. 4a. This shows an average value of 158 Pa m⁻¹. Particularly, large areas with high stress gradients are located below the Hochvernagtspitze in a steep and narrow area. Low stress gradients are largely present below the Gepatschjoch.

Fig. 4. (a) Driving stress gradients in Pa m⁻¹ on Vernagtferner with the areas of high crevasse index from the reference dataset outlined. (b) Identification of crevassed areas using the stress gradient approach with a threshold value of 236 Pa m⁻¹.

Within the crevassed areas, driving stress gradients are significantly higher (average value 224 Pa  m⁻¹) than on the total glacier (average value 158 Pa m⁻¹). This corresponds with our assumptions. The linkage between spatial driving stress gradients and the occurrence of crevasses allows us to investigate the potential of defining crevasse prone areas from the stress distribution.

For this purpose, a threshold value for the driving stress gradient is required, above which an area is considered as crevassed. In the literature, no threshold value has been proposed so far. Figure 5 indicates the percentage of glacier and crevassed areas that exceed each selected threshold. The blue curve visualizes the distribution of the glacier area in relation to the driving stress gradients. The purple curve represents the distribution of the previously defined crevassed areas in relation to the driving stress gradients. The threshold above which an ice surface is classified as crevassed should be a value for which the difference between the glacier area (blue curve) and the defined crevassed area (purple curve) is large. This criterion ensures that many crevassed areas are classified in relation to the glacier area, thus minimizing misclassification. The greatest vertical distance between the two lines is mathematically calculated at a threshold value of 236 Pa/m. Figure 4b shows the four possible results false positives (FP), true positives (TP), false negatives (FN) and true negatives (TN) with this threshold value. At this value, 55% of the crevassed areas are classified correctly and are therefore true positive. However, 25% of the overall glacier is classified as a crevassed area, although the real occurrence is only 6%. This means that 22% (25%(positives) −3% (true positives)) of the glacier were misclassified as crevassed areas and are false positives. The overestimation of crevassed areas is 3.7 times the actual crevassed area.

Fig. 5. Selected threshold value for the spatial change in driving stress.

4.2 Curvature

A strong vertical curvature along the steepest slope supports extending flow conditions and should therefore be relevant for crevasse formation. The values calculated for the 2^nd slope line derivative using the ArcGIS (version 10.8.1) tool for vertical curvature are shown in Fig. 6a. Negative values describe convex and positive values concave areas. Convex values are dominant in the areas of high crevasse index from the reference dataset.

Fig. 6. (a) 2^nd slope line derivative with the areas of high crevasse index from the reference dataset outlined. Negative values describe a convex, positive values a concave surface curvature. (b) Identification of crevassed areas using the curvature approach with a threshold value of −0.095.

Class limits must be set to divide the glacier into ‘concave’, ‘convex’ and ‘flat’. In particular, determining the threshold between convex and flat is important for comparing the convex areas of the glacier with the areas of a high crevasse index from the reference dataset. For this purpose, a graph was plotted to indicate the percentage of the areas of high crevasse index in the reference dataset (purple curve) that is defined as convex using different threshold values (Fig. 7). The graph also shows the percentage of the complete glacier (blue curve) that is detected as crevassed for different thresholds and in order to indicate misclassification. When determining the best threshold, a maximum of the areas with a high crevasse index in the reference dataset should be detected, while as little of the entire glacier as possible should be considered convex to ensure low misclassification. To meet these conditions, the largest vertical difference between the percentages of detected crevassed areas and the areas of the entire glacier which are classified as convex is calculated. The resulting threshold value of −0.095 correctly detects 48% of the crevassed areas as TP and classifies 17% of the total glacier. Considering the correctly classified areas, the misclassification as false positives is 14% (Table 1). Figure 6b visualizes the four possible results false positives (FP), true positives (TP), false negatives (FN) and true negatives (TN) for the threshold value −0.095.

Fig. 7. Selected threshold value for the curvature.

Table 1. Compilation of the results with the respective approaches, with CI meaning crevasse index.

4.3 Combination

A correlation analysis was performed to evaluate a possible relation between driving stress and curvature. The results in Fig. 8a show no strong correlation, but they display a certain pattern: crevassed areas primarily have a negative curvature and a positive spatial change in driving stress and, thus, the purple points are mostly located in the bottom right quadrant.

Fig. 8. (a) Correlational analysis of curvature and driving stress gradients for the glacier and the crevassed areas. (b) Correlation analysis with a best fit line of the crevassed areas in the bottom right quadrant and a distance due to spreading of 0.1.e.

A best fit line can be calculated by the crevassed areas in the bottom right quadrant, to investigate the relation between the curvature and the spatial change in driving stress for the crevassed areas. A possible spreading of 0.1 is assumed in order to detect as many crevassed areas as possible with a minimum overestimation.

All points inside this zone are classified, as shown in Fig. 8b. The misclassification is clearly dominant. This approach detects 33% of the areas with a high crevasse index, with 17% of the areas with a less crevasse index detected, resulting in an overestimation of 2.8 times the actual crevassed area (Table I).

5. Discussion and conclusion

The results of the analysis indicate that crevasses appear mainly in extensional flow or in convex areas, since the crevassed areas are primarily located in the bottom right quadrant of Fig. 8a. Table I summarizes the results with the different approaches. For this purpose, the statistical values for recall and precision were calculated in percent. The recall indicates the percentage of crevassed areas that were found correctly, i.e. the percentage of TP out of the manually detected crevassed areas (recall = TP/(TP+FN)). The precision defines the percentage of the reported crevassed areas that are indeed crevassed, i.e. the percentage of TP out of all values found with the algorithm (precision=TP/(TP+FP). In addition, the misclassification, i.e. the percentage of false positives in relation to the complete glacier, is shown. In order to provide a good comparison, the ratio of the percentage of false positives to the recall as well as the overestimation of crevassed areas is shown. The ratio of the percentage of false positive (misclassification) to the recall (detected crevassed areas) shows that the curvature and a combination of driving stress and curvature perform best. However, the combination of driving stress and curvature results in a lower recall, with a value of only 29%. As a result, the parameter curvature is preferable due to its significantly higher recall compared to the combination. However, it must be noted that for the parameter curvature, the vertical curvature was calculated along the steepest slope and thus transverse crevasses are more likely to be found. In the study area, this is clearly the predominant crevasse type.

However, when considering misclassifications, it should be noted that these values are based on the areas indicated as crevassed areas due to a high crevasse index in the reference dataset. Therefore, it is not necessarily a misclassification in locality, especially if this misclassification is close to the areas detected as crevassed areas. In this case, the crevassing is less than in the areas with a high crevasse index from the reference dataset. Furthermore, as pointed out by Nath and Vaughan (2003), there may be crevasses below the surface that were not detected in the optical reference dataset. This is why the result will actually be slightly better than shown in the analysis, since many misclassified areas are close to the areas indicated as crevassed and are, therefore, probably still actually crevassed, but with a smaller fraction of crevasses. It must also be taken into account that crevasses are interacting with surrounding crevasses, which might change the stress field (Sassolas and others, 1996). Since we investigate crevassed areas, we neglect this fact.

The results point to the probability that it might be easier to identify non-crevassed areas with less misclassification compared to identifying crevassed areas. This may help in applications for alpinists to define safe routes across glaciers, based on the derivation of a safety index, noting that false-negatives are most dangerous for field parties, as they falsely suggest that crevassed areas are not crevassed.

We are aware that our research has three limitations. First, it can only be applied to alpine glaciers with low advection, where crevasse patterns are relatively stable and represent the local strain. Otherwise crevasse pattern may indicate deformation history rather than actual strain (van der Veen, 1999; Lawson and others, 2000; Jennings and Hambrey, 2021). Second, the ice thickness must be known. Third, there is a significant amount of misclassification. These limitations hinder the potential of elevation and ice thickness data for deriving crevassed areas using only a simple parameter-based approach. Further data collection from other glaciers is required to determine how well the change in driving stress and the curvature represents crevassed areas in other environments. In addition, a large data set could be subjected to statistical analysis to evaluate the results.