Remote Sensing Based Land Cover Map of Watersheds in the Swedish Arctics: Study of Spatial and Temporal Variabilities of Land Cover

1. Introduction

When it comes to hydrological and biogeochemical fluxes on continental surfaces, the nature of the land cover, including both vegetation covers, bare rock outcrops and surface water bodies is of major importance [1,2,3]. It is especially true in the Arctic, where permafrost conditions may strongly control the present ecotypes and their distributions [4,5,6,7], while vegetation cover variability may in turn significantly impact thermo-hydrological conditions through evapotranspiration for instance [8,9]. Thus permafrost modeling requires detailed knowledge of the land cover distribution.

The vast extension and the remoteness of the Arctic regions make the establishment of field survey – based land cover maps difficult. Moreover, high resolutions and open data maps are needed for many applications [10]. Thus there is a growing interest in airborne [11,12] and remote sensing [13,14,15] capable of producing high-resolution vegetation maps in the Arctics. These regions are experiencing intensive climate change [16]. Permafrost thawing results in methane emissions [17] which contributes to the greenhouse effect. These modifications induce changes in ecotypes [18] that are visible at the landscape level. Thus there is a need for not only high spatial resolution maps, but also for high temporal resolution survey. In order to produce regularly updated land cover map for large areas, the use of remote sensing data from long term satellite missions combined with in-situ information is required [14].

Here we present a workflow for making high resolution vegetation map using only open data and open source softwares along with dedicated field data. The workflow is applied to two watersheds in the Swedish Arctic, for parts of which previous vegetation maps at lower resolutions and/or in past climatic conditions where already available [11,19,20]. The obtained map is used for studying link between topography and vegetation distribution, and also the temporal evolution of the vegetation cover during a 14 years period (2008-2022).

2. Materials and Methods

2.1. Field sites

Two watersheds close to Abisko Scientific Research Station are considered (Figure 1). The first one from West to East is Miellajokka, a sub-alpine catchment which include the iconic mounts of Tjuonavagge (Lapporten). This 51.5 km² catchment presents altitudes ranging from 383 to 1731 m above see level [21]. The most eastern watershed is Stordalen, a 16 km² catchment with a lake-rich, peat-rich Northern part, and a sub-alpine Southern part, with elevation between 350 and 770 m above see level [22,23,24,25]. In Stordalen vegetation maps of the Northern, low elevation part has been already produced with airborne data of 2000 [19]. Later on, another vegetation map for the whole watershed has been produced with airborne data of 2008 [20]. Both Stordalen and Miellajokka are encompassed in the area studied by Reese, with a vegetation map established on the basis of 2010 satellite images, using also data acquired by lidar survey [11].

2.2. Satellite data

Obtaining images in the Arctic zone to study vegetation cover is difficult. These geographical areas are covered with snow for a large part of the year, which prevents any satellite study of the vegetation cover. In addition, frequent clouds hinder the acquisition of optical images. A single Sentinel-2 image acquired on 25 August 2022 is downloaded from https://peps.cnes.fr. Ten bands were selected for land cover classification (Table 1). The image is not corrected for atmospheric effects (Level-1C). The images are stored in the UTM34N reference coordinate system and all calculations are performed in this system to avoid altering the radiometry by re-projection.

On the basis of four (B03, B04, B08, B11) out of the ten acquired channels, four derived indices are calculated: Bright, NDVI, NDWI, NDII (Table 2). The bright index is very sensitive to albedo. It distinguishes between light and dark soils. The NDWI (Normalized Difference Water Index) was used to detect water areas. The NDVI expressed the photosynthesis of the vegetation cover. The use of NDII [26,27] did not improve the results and was not retained for the final classification.

Since vegetation in mountainous areas is related to altitude, the digital terrain model is a very useful data source. ArcticDEM is an NGA-NSF public-private initiative to automatically produce a high-resolution digital surface model of the Arctic using optical stereo imagery. The majority of ArcticDEM data was generated from the panchromatic bands of the WorldView-1, WorldView-2, and WorldView-3 satellites and, for a small percentage of data, from the GeoEye-1 satellite. For this study, ArticDEM Release 7 "mosaic" format files with a spatial resolution of 2 m are downloaded at https://data.pgc.umn.edu/elev/dem/setsm/ArcticDEM/mosaic/v3.0/.

2.3. Field survey

The field mission took place from 21 July 2022 to 24 July 2022 in the Miellajokka and Stordalen watersheds, northern Sweden. A "Samsung Galaxy Tab S6 Lite" tablet is used to perform the field surveys. It supports GPS, GLONASS, Beidou and Galileo navigation systems. The Qfield software is used for data entry in the field. Its compatibility with QGIS simplifies data collection and subsequent analysis [28].

A database including a color composite of Sentinel-2 image channels B08/B04/B03, the Open Street Map data and an empty vector layer intended for the field surveys is prepared in QGIS then transferred to Qfield. The surveys are made in the form of polygons drawn once in the field at the polygon locations with the pen of the tablet. Each polygon is associated with the observed land use class at the considered place, and possibly a photo is taken with the tablet camera. The land cover classes are initially defined from Reese [11]. However, the 270 surveys conducted during our mission only identified 7 classes (Table 3) out of the 12 classes of Reese. The rock class is mainly composed of fresh bed rock outcrops but it may also includes thin organic soil and sediment. The alpine meadow was not encountered enough during the field trip to constitute a class. Likewise, the Mountain birch - meadow type class was only encountered in 6 surveys and is grouped with the Mountain birch - moss type class to form a single Mountain birch class. Areas in snow that were very poorly represented, during the July mission, are not included. Grass heaths were not encountered. A new class, 'Human infrastructure', is added, representing mainly the road and the railway passing through the mapped area. Shadows in the steep areas to the south of the study area hinder recognition of the landscape they cover. To avoid confusion, especially with water, a shadow class is created, making the number of considered classes equaled to 9.

2.4. Classification

The Random Forest method [29] uses decision trees and random draws without replacement of samples and variables to classify the Sentinel-2 image. Within each class, 30% of the polygons are randomly drawn and reserved for classification quality assessment. The classification is trained with the remaining 70%. GRASS software is used for the calculations [30]. The extension r.learn.ml2 interfaces with the Scikit-learn library written in python to perform Random Forest classification.

3. Results

3.1. Vegetation map in current climatic conditions

The classified image (Figure 2) shows patterns consistent with the knowledge of the terrain. The transport infrastructures are described with precision in their continuity. The lake Torneträsk at the North is homogeneously identified.

The statistics confirm the quality of the classification. The percentage of pixels correctly classified by Random Forest is 83%. The confusion matrix (Table 4) provides an analysis of the accuracy of the classification used for building our map at the class level. The columns show the field surveys and the rows show the classification results. For instance, the number 117 corresponds at the cross of the “mountain birch” line and the “mesic heath” column means that 117 pixels that have been classified as “mountain birch” belong to “mesic heath” according to the field survey. The largely dominant weight of the diagonal terms in this confusion matrix demonstrates the quality of the classification. The confusion between “Dry heath” and “Mesic heath” is understandable because these two formations are differentiated primarily by canopy height, a feature not accessible from the images used in our study. Likewise, the confusion between “Alpine willow” and “Wetland” is due to the difficulty of recognizing spaces occupied by a few willow plants. Besides, with such a pixel classification approach, places that are temporarily flooded at the moment of the satellite image acquisition are difficult to distinguished from true wetlands, i.e. places that are under water almost all along the active season. This could lead to overestimate the wetland area, including into it places with other vegetation types like meadow temporarily flooded by ground water discharge or snowmelt water. Another important point is the detection of temporary high elevation open water bodies in several places around Lapporten lake, according to both this classification and the two indices NDII, NDVI values of the pixels. These ones may be generated by late snowmelt in the highest places of the landscape. Finally, the confusion between “Dry heath”, “Mesic heath” and “Mountain birch” are classes that can be contiguous and even associated in some places. It describes mixed spaces where several classes coexist, i.e. ecotone between these classes.

The “Mountain birch” class is the most represented (Table 5). It occupies more than a quarter of the space. The rock class and the wetlands (“Wetland” and “Alpine willow”) share 40% of the spaces. The other classes are much smaller.

3.2. Influence of altitude on land cover

Land cover appears to be strongly conditioned by altitude. Altitudinal zonation of land cover is encountered in various arctic contexts [31,32,33], including in the Abisko region [34,35]. In Miellajokka catchment, the seasonal variability of the hydrogeochemistry of the stream indicates a strong altitudinal control on hydrological processes, especially during the spring freshet [36], while hydrological conditions strongly interact with vegetation and CO₂ fluxes [37]. Table 6, constructed by cross-referencing the land cover map with the ArticDEM, illustrates this phenomenon in the studied watersheds.

The subalpine level at altitudes below 600 m is mainly occupied by the “Mountain birch” class. The Alpine stage includes the “Dry heath”, “Mesic heath”, “Wetland” and “Alpine willow” formations. It is divided into two sub-stages: Between 600 and 800 m of altitude, the “Mountain birch” class is still very present. Above 800 m, this class gives way to the rock. The nival stage is composed only of rock, probably because it depends on harsher life conditions and more intense erosive processes at higher elevations.

3.3. Comparison between past and present vegetation maps

Three maps of parts of our study area have been produced by different authors. A map was constructed from aerial images of 8 August 1970 and 29 July 2000 [19], but the area covered is too small to allow comparison with our data. Another map of the Stordalen watershed was obtained from images obtained from a helicopter flight on 1 August 2008 [23]. The most recent map of the Miellajokka watershed was produced from SPOT5 images of 28 July 2011 and laser data acquired under leaf-on conditions from two scanning dates (20 August 2010 and 9 September 2010) [11,38]. As the semantics of these maps are not identical to ours, an analysis of the variable typologies is necessary prior to the study of the landscape evolution. In order to overlay the maps and then calculate statistics, the Lundin and Reese maps extracted from the publications were georeferenced from control points identified in the landscape.

Comparison with the map of Reese

Reese [11,38] produced a land cover map with a larger number of classes. In order to compare the Reese map with our data, some classes of the Reese map are combined. “Snow ice” and “Snow bed” classes are grouped together. “Dry heath”, “Extremely dry heath”, “Grass heath” are also grouped together. The “Human infrastructure” class is not considered because it does not exist in Reese's work. "Alpine meadow” and “Tall Alpine meadow" were not confirmed by ground observations during our field trip.

The spatial distribution of the classes is slightly different (Table 7 and Table 8). The “Rock” class accounts for 36% for our classification and only 14% for the Reese map.

The change matrix shown in Table 9 encompasses the differences in semantics of classes and the changes in the landscape between the dates of the two maps (i.e.: 2014 and 2022), and may be also discrepancies due to the use of different methodologies. Table 9 and Figure 3 show that some “Rock” areas in the middle our map are covered by “Grass heath”, “Dry heath” and “Alpine willow” on the Reese map. The forest is also growing slightly to the south in sparse patches.

On the other hand, “Alpine willow” is also more represented on Reese map without any conclusion being drawn because this class is misclassified by Reese [11]: the confusion matrix indicates 20 of 44 pixels misclassified. In addition, “Alpine meadow” which we have not taken into account, is identified as “Wetland” on our map (Figure 3), maybe because most of “Alpine meadow” places were temporarily flooded at the time of observation(see also Section 3.1). This class “Wetland” exists on the map published by Borgelt [38] but it is not present in the confusion matrix published by Reese [11]. The other classes do not show significant difference in in their proportion and spatial distribution.

Despite the differences between the two maps, Figure 3 shows a great consistency in the distribution of space. Three elements can explain these differences. 1) The landscape is natural. There are no parcels to structure it. Between areas occupied by two vegetation classes there is often a transition zone which is difficult to assign to a class. 2) The class definition of Reese [11] takes into account the height of the stratum using metrics derived from laser acquisitions, a technology we did not employ. 3) The landscape has evolved between 2010 and 2022. For example, some areas of rock have become vegetated.

Comparison with the map of Lundin

The definition of classes in Lundin [23] is very different from the one presented in this work. It is therefore possible to compare these two maps. Grouping “Dry heath”, “Mesic heath”, “Alpine willow” allows them to be compared areas to the “Alpine tundra” class of Lundin. Similarly, our “Wetland” class is compared to the “Peatland” class of Lundin. The classes “Human infrastructure” and “Non vegetated” correspond to road, railway and building. The latter class is much more represented on the Lundin map (Table 10 and Table 11). The differences are related to a larger road and railway footprint, which does not affect the landscape dynamics. This observation shows the satisfactory overlay of the two maps.

The change matrix (Table 12) shows that the "Rock" class is also more represented on the map of Lundin, which covers twice as much space. The areas of this class that are not classified as rock for our map are located in the south of the map (Figure 4). They are contiguous to the “Rock” areas of our map. It is also possible that there has been a forest expansion between 2008 and 2002 in this area. Indeed, the forest has grown between 2008 and 2022. It has gained some space in all land cover categories. Furthermore, some areas of “Peatland” appear to be transformed into “Wetland” but Lundin [23] indicate "Peatlands were subdivided into wet areas (fen) and dry areas (bog) proportionally to what was found by Malmer et al. (2005)". It is therefore not possible to draw a conclusion from this observation. In summary, it is possible to compare our map with the map of Lundin after a semantic analysis of the categories. The main differences between the two maps concerns the south of the Stordalen watershed where Alpine toundra and forest are intermixed. The forest seems to have taken over areas previously occupied by tundra.

4. Discussion

The methodology implemented requires little sampling effort to quickly obtain a land cover map. This operation is feasible every year to monitor the evolution of land use. These results are particularly important in a region subject to climate change which is currently undergoing major upheaval.

The results of the Sentinel-2 image classification show a structuring of the landscape as a function of elevation, with different vegetation level along with altitude [39]. In the Miellajokka and Stordalen catchments, three levels are present. The subalpine level (< 600 m) is mainly occupied by “Mountain birch” formations. The alpine level [600 m , 1100 m] is characterized by “Dry heath”, “Mesic heath” and “Alpine willow”. This level could be split into two sub-levels according to the respective abundance of “Moutain birch” at lower altitude and “Rock” at upper altitude. The upper level, the nival level (> 1100 m), is only composed of rock.

Comparison between past and present vegetation maps is made difficult by the lack of a common typology field survey protocols. Nevertheless, It was possible to see a change in the landscape between 2008 and 2022. The analysis of the map of Lundin [23] and the one presented in this work would show an extension of the forest on the tundra towards the south (i.e., toward higher elevations) during the 2008-2022 period. This finding is in agreement with the map of Reese [11] which shows also a slight expansion of forest between 2011 and 2022 in the same direction. This trend could be a consequence of the on-going climate warming.

The proximity of the Abisko observation station makes the Meillajokka and Stordalen watersheds a privileged study area for the evolution of landscapes in the Arctic zone, in particular the thawing of the permafrost at high altitudes which is attested. This monitoring requires annual surveys according to a unified protocol in terms of sampling effort, definition of classes and method of recording the surveys to monitor the evolution of this region subject to climate change.