Monitoring Urban Green Infrastructure Changes and Impact on Habitat Connectivity Using High-Resolution Satellite Data

: In recent decades, the City of Stockholm, Sweden, has grown substantially and is now the largest city in Scandinavia. Recent urban growth is placing pressure on green areas within and around the city. In order to protect biodiversity and ecosystem services, green infrastructure is part of Stockholm municipal planning. This research quantiﬁes land-cover change in the City of Stockholm between 2003 and 2018 and examines what impact urban growth has had on its green infrastructure. Two 2018 WorldView-2 images and three 2003 QuickBird-2 images were used to produce classiﬁcations of 11 land-cover types using object-based image analysis and a support vector machine algorithm with spectral, geometric and texture features. The classiﬁcation accuracies reached over 90% and the results were used in calculations and comparisons to determine the impact of urban growth in Stockholm between 2003 and 2018, including the generation of land-cover change statistics in relation to administrative boundaries and green infrastructure. For components of the green infrastructure, i.e., habitat networks for selected sensitive species, habitat network analysis for the European crested tit ( Lophophanes cristatus ) and common toad ( Bufo bufo ) was performed. Between 2003 and 2018, urban areas increased by approximately 4% while green areas decreased by 2% in comparison with their 2003 areal amounts. The most signiﬁcant urban growth occurred through expansion of the transport network, paved surfaces and construction areas which increased by 12%, mainly at the expense of grassland and coniferous forest. Examination of urban growth within the green infrastructure indicated that most land area was lost in dispersal zones (28 ha) while the highest percent change was within habitat for species of conservation concern (14%). The habitat network analysis revealed that overall connectivity decreased slightly through patch fragmentation and areal loss mainly caused by road expansion on the outskirts of the city. The habitat network analysis also revealed which habitat areas are well-connected and which are most vulnerable. These results can assist policymakers and planners in their e ﬀ orts to ensure sustainable urban development including sustaining biodiversity in the City of Stockholm.


Introduction
Over the past two decades, the City of Stockholm, Sweden, has grown substantially and is now the largest city in Scandinavia. Its population has risen by over 28% since the year 2000. Urbanization due to population increase is placing pressure on the green areas in and around Stockholm, which are Introduction [5,63], often formed as "green wedges", i.e., large swaths of natural and seminatural vegetation that lead from the surrounding rural areas in towards the more densely built-up urban city center. The wedges often contain nature reserves or national parks, such as Stockholm's National City Park, and are part of Stockholm regional planning to ensure recreational opportunities, habitat for wildlife and the provision of ecosystem services [63]. Yet, the City's natural areas are under pressure from the region's growing population and the increasing demand for built-up areas. Plans for the creation of 22,000 new dwellings each year up to 2030 are already underway [63]. This study investigates the impact of urban expansion on green areas and green infrastructure in Stockholm City between 2003 and 2018.
Stockholm City is the study area and comprises 14 city districts ( Figure 1). The total analysis area is approximately 216 km 2 including water. Land-cover classes in the city that could be identified from the satellite imagery with the selected classification techniques are buildings, transport/construction sites, urban green space, grass/open fields, coniferous forest, deciduous forest, rock outcrop, barren land/soil, agricultural land, wetlands and water. The transport/construction class includes roads, railway, airport runway and pavement in general, as well as construction and industrial extraction sites. Prime examples of urban green space are city parks, allotment areas or gardens, but the class also includes any otherwise undesignated green area composed of a mix of trees, shrubs and grass in or around urban areas.  Two WorldView-2 (WV-2) images at product level Standard 2A captured on 4 July 2018 and covering slightly more than the geographic area of Stockholm City were acquired for classification. Two scenes of QuickBird-2 (QB-2) images at product level Standard 2A captured on 14 June 2003 (covering approx. 90% of Stockholm City) and 6 September 2005 (covering approx. 10% of Stockholm City) were also acquired. Efforts were made to obtain images from the peak of the vegetation season to take advantage of spectral differences between green and gray areas and in order to avoid the detection of nonpermanent changes caused by seasonal differences. The same four multispectral bands from both sets of imagery were used for segmentation and classification, namely the blue, green, red and near-infrared bands. Only four of the eight multispectral bands available in the WV-2 data were used in order to reduce processing time and to not exceed computer memory constraints as well as to match the spectral input for the segmentation and classification of the 2003 QB-2 imagery.
The City of Stockholm has invested in the development of geographic information about its green infrastructure in order to facilitate the integration of biodiversity aspects into planning and environmental assessment. Stockholm City developed a geographic dataset representing its green infrastructure [64]. This polygon vector dataset was created in 2013 based on the city's geographic ecological information sources with reference years varying between 1995 and 2010, such as the oak database, the city biotope map and its species observation data archive (ArtArken [65]), as well as on representative habitat networks. These habitat networks represent selected sensitive species and were developed by Mörtberg et al. [66,67], based on available geographic data and expert knowledge. Likely impacts on the habitat networks were estimated based on a then-current city planning scenario, but they have not yet been used in comparison with spatiotemporal assessment of landscape changes from urban growth.
Ecologically significant areas in the green infrastructure dataset are divided into three categories, referred to as ecologically significant (1) core areas, (2) habitat for species of conservation concern and (3) dispersal zones. The ecological qualities of the core areas mean that they often host a wide variety of plant and animal species and play an important role in maintaining ecosystem functions. Habitat for species of conservation concern is based on occurrences of red-listed species and their habitat requirements. Dispersal zones are areas important for the movement of species between habitat patches or core areas.

Image Preprocessing
An overview of the study methodology is presented in Figure 2. Processing Standard 2A provides WV-2 and QB-2 images that are radiometrically corrected and orthorectified. Thus, only two preprocessing steps were necessary. First, the WV-2 and QB-2 multispectral bands (spatial resolutions: 2 and 2.4 m, respectively) were pansharpened to 0.5 and 0.6 m, respectively, based on their panchromatic bands. The higher spatial details of the pansharpened images were deemed essential for improved classification of small houses and small urban green spaces such as allotment areas and gardens [68]. The pansharpening was performed using the automatic image fusion algorithm in PCI Geomatica developed by Zhang [69]. Secondly, since two scenes of both types of imagery were required to cover the geographic area of Stockholm City, both pairs of images were mosaicked, also using PCI Geomatica. Once these steps were completed, however, a slight shift (around 5-10 m) due to a residual geolocation error was discovered between the QB-2 and WV-2 images. Given the higher geolocation accuracy of the WV-2 sensor (CE90 < 3.5 m) compared to QB-2's (CE90 23 m) [70,71], we used the WV-2 image as reference and coregistered the images using the automatic coregistration image workflow available in ENVI version 5.3 (Exelis Visual Information Solutions, Boulder, Colorado). The final images were correctly coregistered with subpixel residual error. Remote Sens. 2020, 12, x FOR PEER REVIEW 6 of 30

Segmentation and Classification of Satellite Data
Object-based image analysis (OBIA) has become widely used in recent years as a key step in the classification of high-resolution data thanks to its frequently superior performance than many pixel-based approaches [19,72,73]. OBIA was deemed appropriate for this study given the resolution of the data and the focus of the classifications, which is on nonurban vs. urban land-cover types. eCognition Developer software [74] was used to segment, extract features from and classify the two sets of images. Scale parameters of 100 for 2003 and 125 for 2018 were used for segmentation. The homogeneity criteria of shape and compactness were set at 0.1 and 0.5, respectively, for both years. Parameter selection was based on trials and yielded the most appropriate objects for classification by virtue of their close representation of discrete areas of various land/cover types and structures in the images.
Training samples were gathered for 17 classes initially: 6 classes for buildings based on varying reflectance properties of building materials (for example, copper roofing generally had greater reflectance values than concrete roofing), 2 for transport lines and construction sites, 1 for urban green space, 2 for open grassland and grass sports fields, 1 for coniferous forest, 1 for deciduous forest, 1 for bare rock, 1 for barren land/bare soil, 1 for water and 1 for shadow. Following classification, spatial relationships between classes such as relative border thresholding were used to correct object misclassifications. In this way, many isolated misclassifications were corrected within the coniferous forest, deciduous forest and water classes in particular. Objects identified as shadow were assigned to the class with which they shared the largest relative border. Roads misclassified as buildings were corrected using a length threshold.
Classification of a wetland category was tested but it was heavily confused with both the coniferous forest and rock outcrop classes. Wetlands cover less than 1% of Stockholm City. Given their limited coverage and significant confusion with other classes, wetland areas were masked in

Segmentation and Classification of Satellite Data
Object-based image analysis (OBIA) has become widely used in recent years as a key step in the classification of high-resolution data thanks to its frequently superior performance than many pixel-based approaches [19,72,73]. OBIA was deemed appropriate for this study given the resolution of the data and the focus of the classifications, which is on nonurban vs. urban land-cover types. eCognition Developer software [74] was used to segment, extract features from and classify the two sets of images. Scale parameters of 100 for 2003 and 125 for 2018 were used for segmentation. The homogeneity criteria of shape and compactness were set at 0.1 and 0.5, respectively, for both years. Parameter selection was based on trials and yielded the most appropriate objects for classification by virtue of their close representation of discrete areas of various land/cover types and structures in the images.
Training samples were gathered for 17 classes initially: 6 classes for buildings based on varying reflectance properties of building materials (for example, copper roofing generally had greater reflectance values than concrete roofing), 2 for transport lines and construction sites, 1 for urban green space, 2 for open grassland and grass sports fields, 1 for coniferous forest, 1 for deciduous forest, 1 for bare rock, 1 for barren land/bare soil, 1 for water and 1 for shadow. Following classification, spatial relationships between classes such as relative border thresholding were used to correct object misclassifications. In this way, many isolated misclassifications were corrected within the coniferous forest, deciduous forest and water classes in particular. Objects identified as shadow were assigned to the class with which they shared the largest relative border. Roads misclassified as buildings were corrected using a length threshold.
Classification of a wetland category was tested but it was heavily confused with both the coniferous forest and rock outcrop classes. Wetlands cover less than 1% of Stockholm City. Given their limited coverage and significant confusion with other classes, wetland areas were masked in based on the Remote Sens. 2020, 12, 3072 7 of 29 corresponding class in a 10 m national land-cover dataset [75] and checked manually post-classification. Agricultural fields were also incorporated as a separate class in this way, using a mask based on the appropriate category of the national land-cover dataset.
Various object features and classification algorithms were tested and compared. A successful approach used in previous research [76] with slightly lower resolution imagery was adopted here and yielded good results. This involved a support vector machine (SVM) classification algorithm trained with statistics from one textural, two spectral and three shape data features. Spectral band mean and standard deviation were used as well as Roundness, Asymmetry and Rectangular fit shape features and the Gray-Level Difference Vector (GLDV) Entropy (all directions) texture feature contributed to improved classification results.
To further improve the classification results and their comparability, corrective rules were used under urban and forest masks to reclassify noted problematic areas. A low-density built-up area mask from a previous 2005 classification was manually verified/edited and used to standardize the classifications for building and transport areas that had not changed. Unchanged major transport lines from 2003 were compared to those in the 2018 classification and added where necessary. A forest mask based on the relevant classes from the national land-cover dataset [75] was applied to improve deciduous forest classification where shadow was confused with coniferous forest in the 2018 classification. Green roofing added during the study period was controlled for through an overlay check of conversion from urban to green structure classes and misclassifications were manually corrected. Once post-classification corrections were completed, a minimum of 500 validation sample vector points for each land-cover class were used for accuracy assessment. The points were manually recorded for both years with the satellite imagery, historical Google Earth images and local knowledge of the study area as comparing sources.

Landscape Change Analysis
Land-Cover Change According to Administrative Boundaries and Within the Green Infrastructure Based on the 2003 and 2018 classifications, the land area percent change in urban areas (comprised of the buildings and transport/construction classes) and in green areas (comprised of deciduous/coniferous forest, urban green space, grasslands, agricultural land, barren land, wetlands and bare rock) per city district was calculated. This information reveals where and relatively how much change occurred in different areas within the city. Urban land-cover change in relation to the three types of ecologically significant areas of the green infrastructure as identified by the city was also investigated. The increase in the amount and percentage of urban overlap in these areas was quantified and the change shown graphically.

Change in Habitat Network Connectivity
Habitat connectivity plays a significant role in biodiversity conservation due to its effect on metapopulation dynamics [45], even if habitat amount has been identified as the major contributing factor to species persistence and abundance [43,44]. Habitat networks for selected species sensitive to urbanization were developed by Mörtberg et al. [66,67], including the European crested tit (Lophophanes cristatus) and the common toad (Bufo bufo).
The crested tit is strongly tied to mature and old pine and spruce forests where it requires access to deadwood in standing trees for nesting. It is also sensitive to fragmentation since it avoids open areas and prefers relatively large (10-25 ha), well-connected coniferous forest patches [77][78][79][80]. Information on habitat requirements, gathered by a literature review and interviews with experts by Mörtberg et al. [67], provided the basis for assumptions, criteria and parameters used for the crested tit habitat network analysis in the current study (Appendix A).
Patches for the crested tit habitat analysis were derived by extracting all mature and old coniferous forest from the classifications under a mask of the crested tit habitat geographic dataset created by Mörtberg et al. [67]. All patches with an area of 2 ha or greater were selected initially for analysis since this is below the assumed 3 ha minimum requirement for nesting areas but is still assumed large enough to be used for foraging [67]. A threshold for patch area was necessary for two reasons: to limit the number of patches in order to reduce the required computational resources so that the habitat network analysis could be performed, and because, upon comparison, it was found that larger patches tended to be characterized by higher habitat quality (based on age and type of coniferous forest), one of the main attributes of the crested tit habitat.
Patches that crossed the city boundary were kept in their entirety so as not to unduly influence the network analysis through truncated patches. The preliminary habitat patch selections for 2003 and 2018 were then cross-checked and modified if necessary to ensure that approximately the same patches were included for both years and that none were eliminated because of minor misclassifications or slightly different delineations. Patches with areas down to 1 ha were added where necessary. The resulting habitat networks for both years contain approximately 230 patches each with an average patch size of about 9 ha.
Six of the 13 amphibian species present in Sweden can be found in Stockholm County and Mörtberg et al. [66] selected the common toad as a study object for the development of a habitat network tool. Like other amphibians there, it is dependent on both land and water during its lifecycle and reproduces in water. While the toad is not as threatened by fish since its eggs are poisonous, can travel relatively long distances and is somewhat of a generalist, it is nonetheless sensitive to small barriers and traffic. As highlighted by Mörtberg et al. [66], this makes it suitable as a representative species in landscape ecological assessments of urbanization where habitat fragmentation is addressed and where there is less focus on local habitat variations. While toad populations in Stockholm City have become stable thanks to recent local initiatives, they are nonetheless still under pressure from urban expansion [81,82]. Again, information on dispersal and habitat requirements, as gathered by a literature review and interviews with experts by Mörtberg et al. [66], provided the basis for assumptions, criteria and parameters used for the common toad habitat network analysis in the current study (Appendix A).
Mörtberg et al. [66] created a geographic data layer of common toad reproductive areas in Stockholm City based on proximity to water and appropriate land-cover types. This dataset was compared to urban areas (transport network and buildings) in the 2003 and 2018 classifications and those reproductive locations that overlapped more than 50% were removed. Patches within the City of Stockholm and up to 300 m from its boundaries with an area of 500 m 2 or greater were selected. These areal and spatial limitations were necessary in order to reduce the required computational resources so that the habitat network analysis could be performed. This yielded a 2003 reproductive area dataset of 838 patches with an average patch size of 0.04 ha and a 2018 dataset of 811 patches also with an average patch size of 0.04 ha.
The next step was to construct a friction surface representing in numeric form the resistance or degree of difficulty that a toad might experience while moving through the Stockholm City landscape. Two surfaces (for 2003 and 2018) were created based on the land-cover classifications where each land-cover category was associated with a numeric cost. The cost was established based on the dispersal profile for the common toad outlined in Mörtberg et al. [66] who had also created an example friction surface in 2006. The costs for each land-cover type present in the high-resolution classifications as well as from the study by Mörtberg et al. [66] are listed in Table 1.
Based on the friction surfaces, least-cost path networks connecting the reproductive areas were created for both 2003 and 2018. Mörtberg et al. [66] noted that dispersal or movement distances for the common toad often do not exceed 2 km across optimal (low-cost) land-cover. For this reason, cost paths longer than 2 km were removed from the networks. Water 0-150 m from shore friction gradient that increases exponentially with increased distance, from 1 to 10000 Conefor 2.6 is a freely available software package that can calculate a range of network connectivity indices [83] and is here used to calculate the probability of connectivity (PC) metric, using the following function: where a i and a j correspond to habitat attributes of patches i and j, of which the area is used in this study. A L is the maximum landscape attribute and therefore the same as the total landscape area in this case. p * ij is the maximum product probability of dispersal and equals the product of probabilities of traveling through adjacent habitat patches to get to a destination patch. It is used in index calculation if larger than the probability of direct dispersal (p ij ).
The contribution of landscape elements (patches and links) to overall habitat connectivity and availability can be calculated from the percentage of variation in PC (dPC k ) caused by the removal of each single element from the landscape [46,52]: where ∆PC k is the difference in the value of PC when element k is present compared to when it has been removed from the network. Following removal, any maximum product probability paths that previously went via k have to be rerouted through other elements available in the network. Saura and Rubio [53] demonstrated how dPC k can be partitioned into three fractions which represent various ways in which a patch or node can contribute to habitat availability and connectivity in a network: dPCintra k is a measure of internal connectivity and corresponds to a i x a j (Equation (1)) when i = j = k and p * ij = 1. It measures the habitat availability in the patches themselves independent of other elements in the network. dPC f lux k is similar to other measures of flux [46,84,85] but considers maximum product probability paths (p * ij ) instead of probabilities of direct dispersal (p ij ). It corresponds to the case where i = k or j = k, but i j. dPCconnector k depends solely on the topological position of a connecting element in the landscape network and is a measure of a patch or link k's contribution to connectivity between other habitat patches. It is the partial sum of p * ij a i a j (Equation (1)) where i k, j k and k is part of the maximum probability path that lies between them. For a more detailed explanation of the PC metric and its component fractions, see Saura and Rubio [53].
For the crested tit habitat network analysis, patch nodes and the distances between them (measured from patch perimeter) were generated using the ArcGIS extension Conefor Inputs [86] and are used as input files for the calculation of PC and its component fractions. Dispersal or movement distances and their probabilities were selected based on species profile information reported by Mörtberg et al. [67]. Small birds can move easily through a landscape thanks to their ability to fly but given the assumed reluctance of the crested tit to cross open spaces, a median dispersal distance of 200 m was chosen corresponding to probability 0.5. An assumed maximum connection distance of 2000 m was used. Connectivity indices were calculated based on the nodes and distances derived from the coniferous forest habitat datasets extracted from the 2003 and 2018 classifications.
For the common toad habitat network analysis, patch nodes and probabilities of connections between them were used as input files for the calculation of PC and its component fractions. The least-cost path networks provided cost information for each path and in order to convert these to probabilities of usage, the lowest cost available in the network (8.485 in both 2003 and 2018) was divided by the cost of the path in question. This yielded internode path probabilities between 0 and 1. Changes in PC as well as in its component fractions were compared between years for both habitat networks.
An index useful for comparing temporal changes to landscape connectivity can be derived directly from PC. This is the equivalent connected area (ECA) which is defined as the size of one habitat patch that would provide the same connectivity probability value as the habitat pattern in the landscape under examination [54]. It is calculated by taking the square root of the numerator in the PC index (Equation (1)) and has area units. Relative change in ECA (dECA) can also be evaluated and was calculated for this study in the following way: dECA is also of interest since it can be directly compared to dA, which is the change in the total amount of habitat area calculated in the same way. ECA, dA and dECA for the landscapes in the years evaluated are therefore calculated. Table 2 shows the accuracy scores for the segmentation-based SVM classifications after post-classification corrections. A 93% overall accuracy (kappa: 0.92) was achieved for 2003 and over 91% overall accuracy (kappa: 0.90) for 2018. The producer's and user's accuracies show more specifically how well each land-cover type was classified. The accuracies for all classes were above 80% and thus satisfactory. There was some confusion between the building and transport/construction classes and between barren land and transport/construction, which can be expected given how spectrally similar these classes usually are. The coniferous forest class in the 2018 classifications is slightly overestimated at the expense of deciduous forest and urban green space, due in part to longer shadows in the 2018 imagery which were at times classified as coniferous forest within those land-cover classes.   Figure 3.

Urban Land-Cover Change and Environmental Impact Analysis
In relation to their 2003 class areas, buildings increased by 2.3% or approximately 125 ha and transport/construction areas increased by 12% or roughly 202 ha. The classes most impacted by this expansion were coniferous forest, which decreased by approximately 3% or 101 ha, and grassland/open land which were reduced by 10% or 134 ha. Deciduous forest decreased by about 1% or 34 ha. The largest proportional decrease was for barren land which underwent a 19% decrease. This was largely due to the conversion of dirt sports fields to artificial grass, a type of land-cover which was classified as rock outcrop in the 2018 classification and accounts for a small percentage increase in the rock outcrop class. Wetlands and agricultural fields remained relatively stable with a slight increase for wetlands due to the later summer vegetation conditions in the 2018 imagery and a slight decrease in agricultural areas. Urban green space registered a very slight increase of less than 1%.  Figure 4 displays how "urban" each district is with Kungsholmen, Norrmalm and Södermalm as the clear leaders in terms of urban areas with more than 50% urban land-cover and districts like Skarpnäck and Farsta covered by about 30% or less urban land-cover classes. The ratio of urban land-cover to green areas in the majority of the In relation to their 2003 class areas, buildings increased by 2.3% or approximately 125 ha and transport/construction areas increased by 12% or roughly 202 ha. The classes most impacted by this expansion were coniferous forest, which decreased by approximately 3% or 101 ha, and grassland/open land which were reduced by 10% or 134 ha. Deciduous forest decreased by about 1% or 34 ha. The largest proportional decrease was for barren land which underwent a 19% decrease. This was largely due to the conversion of dirt sports fields to artificial grass, a type of land-cover which was classified as rock outcrop in the 2018 classification and accounts for a small percentage increase in the rock outcrop class. Wetlands and agricultural fields remained relatively stable with a slight increase for wetlands due to the later summer vegetation conditions in the 2018 imagery and a slight decrease in agricultural areas. Urban green space registered a very slight increase of less than 1%.  Figure 4 displays how "urban" each district is with Kungsholmen, Norrmalm and Södermalm as the clear leaders in terms of urban areas with more than 50% urban land-cover and districts like Skarpnäck and Farsta covered by about 30% or less urban land-cover classes. The ratio of urban land-cover to green areas in the majority of the remaining districts is roughly 35/65, which indicates that on the whole Stockholm is still a very green city.  The amount of growth in hectares for urban areas and the percentage loss of green areas are listed per district in Table 3. The values are listed from largest to smallest. This listing makes clear that Södermalm experienced the least amount of change over the 15 years, while Bromma had the most hectares converted to urban land-cover and Kista-Rinkeby had the greatest percentage loss of green areas. In Kista-Rinkeby, a good deal of the green area conversion was due to the expansion of the transportation network. For Bromma, this was also true near its airport, although there were instances of building construction as well.   The amount of growth in hectares for urban areas and the percentage loss of green areas are listed per district in Table 3. The values are listed from largest to smallest. This listing makes clear that Södermalm experienced the least amount of change over the 15 years, while Bromma had the most hectares converted to urban land-cover and Kista-Rinkeby had the greatest percentage loss of green areas. In Kista-Rinkeby, a good deal of the green area conversion was due to the expansion of the transportation network. For Bromma, this was also true near its airport, although there were instances of building construction as well.   The amount of growth in hectares for urban areas and the percentage loss of green areas are listed per district in Table 3. The values are listed from largest to smallest. This listing makes clear that Södermalm experienced the least amount of change over the 15 years, while Bromma had the most hectares converted to urban land-cover and Kista-Rinkeby had the greatest percentage loss of green areas. In Kista-Rinkeby, a good deal of the green area conversion was due to the expansion of the transportation network. For Bromma, this was also true near its airport, although there were instances of building construction as well. In total, there was a 9% increase in urban areas within the green infrastructure over the time period. Table 4 reports the areal and percentage increases in urban land-cover in each ecologically significant category of the green infrastructure, also shown in Figure 6. In total, there was a 9% increase in urban areas within the green infrastructure over the time period. Table 4 reports the areal and percentage increases in urban land-cover in each ecologically significant category of the green infrastructure, also shown in Figure 6.    Figure 7 illustrates the location of "newer" urban areas within the green infrastructure that appeared after 2003. These urban areas are shown in yellow while the different ecologically significant categories of the green infrastructure are highlighted in various shades of blue/purple. Roughly two-thirds of the increase within habitat for species of conservation concern was due to the construction of buildings, in areas such as Norra Sköndal (Farsta district). The other third was due to road construction and an increase in paved surfaces occurring, for example, in Solvalla near Bromma airport. Urban growth within core areas was almost exclusively due to the expansion of the transport network, most notably in the core areas located in the Kista-Rinkeby and Spånga-Tensta districts. Dispersal zones were slightly more impacted by transport network expansion, for example by the construction of new piers in the Östermalm district, than they were by new buildings of which there are examples in Beckomberga (Bromma) and Långbro/Solberga (Älvsjö).
Remote Sens. 2020, 12, x FOR PEER REVIEW 14 of 30 Roughly two-thirds of the increase within habitat for species of conservation concern was due to the construction of buildings, in areas such as Norra Sköndal (Farsta district). The other third was due to road construction and an increase in paved surfaces occurring, for example, in Solvalla near Bromma airport. Urban growth within core areas was almost exclusively due to the expansion of the transport network, most notably in the core areas located in the Kista-Rinkeby and Spånga-Tensta districts. Dispersal zones were slightly more impacted by transport network expansion, for example by the construction of new piers in the Östermalm district, than they were by new buildings of which there are examples in Beckomberga (Bromma) and Långbro/Solberga (Älvsjö).

Changes in Habitat Connectivity and Availability
The habitat network analysis revealed that crested tit habitat connectivity and availability decreased over the 15-year period. Table 5 shows the amount of decrease in the connectivity indices between 2003 and 2018.

Changes in Habitat Connectivity and Availability
The habitat network analysis revealed that crested tit habitat connectivity and availability decreased over the 15-year period. Table 5 shows the amount of decrease in the connectivity indices between 2003 and 2018. The overall probability of connectivity decreased slightly for the habitat network. The decrease in ECA signifies that a single habitat patch (maximally connected) providing the same PC value as the habitat pattern in the landscape would in 2018 only be 508 ha as opposed to 530 ha in 2003. The relationship between change in overall habitat area (dA) and the change in the equivalent connected area (dECA) based on the 2003 and 2018 landscapes is dECA < dA < 0 or specifically -0.04 < -0.03 < 0 which meant that the decrease in connectivity was larger than the decrease in the area of habitat patches.
dPC values for all crested tit habitat patches in the 2003 and 2008 landscapes are illustrated in Figure 8. The darker patches are those that are most important in terms of overall connectivity. It is evident that the patch area has a significant influence on dPC value since the largest patches are also the darkest. This is indeed the case since two of the three component fractions of dPC are determined by or strongly connected to patch area: dPCintra measures internal patch connectivity (area) and dPCflux is a measure of area-weighted dispersal probabilities. Areas of high-quality habitat for the crested tit highlighted in the report by Mörtberg et al. [63] are the same as those with relatively higher dPC values in the figures, with the exception of Årstaskogen (on the waterfront in the northern part of the Enskede-Årsta-Vantör district). The number of habitat patches in each value interval is reported in parentheses and there is a slight shift towards lower dPC values in 2018. The overall probability of connectivity decreased slightly for the habitat network. The decrease in ECA signifies that a single habitat patch (maximally connected) providing the same PC value as the habitat pattern in the landscape would in 2018 only be 508 ha as opposed to 530 ha in 2003. The relationship between change in overall habitat area (dA) and the change in the equivalent connected area (dECA) based on the 2003 and 2018 landscapes is dECA < dA < 0 or specifically -0.04 < -0.03 < 0 which meant that the decrease in connectivity was larger than the decrease in the area of habitat patches.
dPC values for all crested tit habitat patches in the 2003 and 2008 landscapes are illustrated in Figure 8. The darker patches are those that are most important in terms of overall connectivity. It is evident that the patch area has a significant influence on dPC value since the largest patches are also the darkest. This is indeed the case since two of the three component fractions of dPC are determined by or strongly connected to patch area: dPCintra measures internal patch connectivity (area) and dPCflux is a measure of area-weighted dispersal probabilities. Areas of high-quality habitat for the crested tit highlighted in the report by Mörtberg et al. [63] are the same as those with relatively higher dPC values in the figures, with the exception of Årstaskogen (on the waterfront in the northern part of the Enskede-Årsta-Vantör district). The number of habitat patches in each value interval is reported in parentheses and there is a slight shift towards lower dPC values in 2018. Habitat patches and least-cost pathways in the common toad habitat network for 2003 (in yellow) and 2018 (in orange) are shown in Figure 9. The habitat patches are in general small (average size 0.04 ha) but the links between them are more visible and the image makes clear that least-cost paths for toads often run through green spaces or along the outer edge of built-up areas. It also shows how certain parts of the network are isolated, notably in the northeast and southwest of Stockholm City. Very low values for PC were obtained for these networks with no noted change in the index, as displayed in Table 6. This is a limitation of PC acknowledged by Saura et al. [54] which can occur when habitat patches are small compared to the total landscape area. ECA as recommended by Saura et al. [54] can overcome this since it has areal units and will not be smaller Habitat patches and least-cost pathways in the common toad habitat network for 2003 (in yellow) and 2018 (in orange) are shown in Figure 9. The habitat patches are in general small (average size 0.04 ha) but the links between them are more visible and the image makes clear that least-cost paths for toads often run through green spaces or along the outer edge of built-up areas. It also shows how certain parts of the network are isolated, notably in the northeast and southwest of Stockholm City. Very low values for PC were obtained for these networks with no noted change in the index, as displayed in Table 6. This is a limitation of PC acknowledged by Saura et al. [54] which can occur when habitat patches are small compared to the total landscape area. ECA as recommended by Saura et al. [54] can overcome this since it has areal units and will not be smaller than the largest patch in the landscape. From Table 6, we see that ECA decreased slightly signifying a drop in connectivity.
Remote Sens. 2020, 12, x FOR PEER REVIEW 16 of 30 than the largest patch in the landscape. From Table 6, we see that ECA decreased slightly signifying a drop in connectivity.  Yet dA decreased more than dECA, which indicates that the loss of habitat area was more significant between 2003 and 2018 than was the loss of connectivity in the network over the same time period. This is likely due to the fact that whole albeit small patches were lost due to the construction of new urban areas, examples of which are shown in Figure 10. Figure 10a shows an area near the Hansta Nature Reserve where 2003 patches (in green) were replaced by the expansion of transport infrastructure, while Figure 10b shows an area near Skarpnäck where 2003 patches were replaced by housing developments. For the crested tit, habitat patches were often reduced in size or divided but did not disappear altogether. The disappearance of patches and associated links in the common toad habitat meant that the network was reshaped and even to some extent  Yet dA decreased more than dECA, which indicates that the loss of habitat area was more significant between 2003 and 2018 than was the loss of connectivity in the network over the same time period. This is likely due to the fact that whole albeit small patches were lost due to the construction of new urban areas, examples of which are shown in Figure 10. Figure 10a shows an area near the Hansta Nature Reserve where 2003 patches (in green) were replaced by the expansion of transport infrastructure, while Figure 10b shows an area near Skarpnäck where 2003 patches were replaced by housing developments. For the crested tit, habitat patches were often reduced in size or divided but did not disappear altogether. The disappearance of patches and associated links in the common toad habitat meant that the network was reshaped and even to some extent consolidated, as illustrated in Figures 9 and 10. Therefore most impact registered through the loss of habitat area rather than connectivity.
Remote Sens. 2020, 12, x FOR PEER REVIEW  17 of 30 consolidated, as illustrated in Figures 9 and 10. Therefore most impact registered through the loss of habitat area rather than connectivity. In order to understand more about the nature of the changes in PC and ECA values in both networks and since dPC is strongly influenced by habitat area, it is of interest to break down the metric into its component fractions and so that the value of patches as connectors or stepping stones can be assessed. Table 7 shows the percent contribution of each component of the PC index (dPCintra, dPCflux and dPCconnector) in 2003, 2018 and the interperiod change for both the crested tit and common toad habitat networks. For the crested tit habitat, dPCflux of habitat patches contributes most to connectivity in the landscape (approximately 50%) given the median distance of 200 m, followed by dPCintra (approximately 27%) and then dPCconnector (approximately 23%). By contrast for the common toad habitat network, dPCintra was the largest contributing factor to PC (approximately 93%) meaning that patches contributed most to connectivity through their size, rather than via flux (approximately 7%) or their connector role (0.2%). It is clear that the patches in the common toad habitat play a very small role as connectors or stepping stones in the network and that the loss in connectivity was due to a small decrease in flux. This differs from how connectivity changed in the crested tit habitat where the decrease was mainly due to patches playing less of a connector role in the network, although the flux was also negatively affected to a lesser degree. Values of dPCflux in the crested tit habitat network are worth examining since this fraction is the largest contributor to overall PC and since the impact on probable connectivity was greater than the impact on habitat area as indicated by the change in ECA (Table 5). Slight decreases in value for some patches were noted in the northwest (Hässelby-Vällingby and Kista-Rinkeby districts) and southeast (southern Skarpnäck) of the city which are due to the construction of new roads resulting in fragmentation of patches in and near Hässelby-Vällingby and loss of habitat area for the two others in Kista-Rinkeby and Skarpnäck. Figure 11 shows the dPCflux values per habitat patch for each year in Hässelby-Vällingby and in southern Skarpnäck. These changes influence the value of nearby patches in the case of dPCflux since decreases in habitat area and increases in distance between patches will affect area-weighted (maximum) probabilities of dispersal. Otherwise, there is an observed consistency of values over most of the remaining landscape.
Examining the link structure from the analysis can reveal areas where restorative conservation efforts may have the most benefit. Figure 12 shows connecting links of up to 500 m in distance with dPC of crested tit habitat patches in 2018. Roughly 20% of crested tits would travel 500 m according to the median distance assumed here (200 m corresponding to probability 0.5). As an example, from this image, we can see that Årstaskogen (white arrow), one of the high-value habitat patches highlighted by Mörtberg et al. [67], appears isolated from the rest of the network. Patches that are exposed and at risk of becoming isolated also become apparent (blue arrow). The most important links for connectivity in the network based on testing patch removal (Equation (2)) were found in the areas of patches that were most well-connected, specifically in Flaten Nature Reserve in southern Skarpnäck and in Bromma. Figure 9 also reveals areas of relatively higher connectivity for the common toad, most notably again in and around Flaten Nature Reserve and Grimsta forest between Hässelby-Vällingby and Bromma. Remote Sens. 2020, 12, x FOR PEER REVIEW 19 of 30 Examining the link structure from the analysis can reveal areas where restorative conservation efforts may have the most benefit. Figure 12 shows connecting links of up to 500 m in distance with dPC of crested tit habitat patches in 2018. Roughly 20% of crested tits would travel 500 m according to the median distance assumed here (200 m corresponding to probability 0.5). As an example, from this image, we can see that Årstaskogen (white arrow), one of the high-value habitat patches highlighted by Mörtberg et al. [67], appears isolated from the rest of the network. Patches that are exposed and at risk of becoming isolated also become apparent (blue arrow). The most important links for connectivity in the network based on testing patch removal (Equation (2)) were found in the areas of patches that were most well-connected, specifically in Flaten Nature Reserve in southern Skarpnäck and in Bromma. Figure 9 also reveals areas of relatively higher connectivity for the common toad, most notably again in and around Flaten Nature Reserve and Grimsta forest between Hässelby-Vällingby and Bromma.  It is made up of several separate groupings of connected habitat patches. Therefore, the probability that two crested tits, if placed randomly within the network, would be able to reach one another is very low at the outset in 2003 (0.06%) and decreased very slightly (-0.01%) by 2018. Recommendations for future planning could include a focus on preserving crested tit habitat particularly in areas where connectivity is still good but where links may be exposed or at risk, such as in south-central Stockholm City and in Bromma. For the common toad, more focus may be necessary on preserving smaller habitat areas since many of the larger habitat patches are already part of the national urban park or nature reserves (such as Grimsta, Judarskogen, Kyrksjölöten, etc.). Attention and planning measures could be directed to preserving the patches that enable links between these areas.  Figure 12 also reveals why the PC values are very low in general for the crested tit habitat network. It is made up of several separate groupings of connected habitat patches. Therefore, the probability that two crested tits, if placed randomly within the network, would be able to reach one another is very low at the outset in 2003 (0.06%) and decreased very slightly (-0.01%) by 2018. Recommendations for future planning could include a focus on preserving crested tit habitat particularly in areas where connectivity is still good but where links may be exposed or at risk, such as in south-central Stockholm City and in Bromma. For the common toad, more focus may be necessary on preserving smaller habitat areas since many of the larger habitat patches are already part of the national urban park or nature reserves (such as Grimsta, Judarskogen, Kyrksjölöten, etc.). Attention and planning measures could be directed to preserving the patches that enable links between these areas.

Discussion
Urban growth in Stockholm City between 2003 and 2018 was modest given the increase in its population over the same time period. In relation to their 2003 class areas, buildings increased by 2.3% and transport/construction areas increased by 12% over the study time period. Coniferous forest and grasslands/open fields were the vegetated land-cover types most impacted by this urban growth, with a loss greater than 100 ha each. Expansion of the transport network negatively affected the core areas of the green infrastructure, while construction of buildings had a greater impact on

Discussion
Urban growth in Stockholm City between 2003 and 2018 was modest given the increase in its population over the same time period. In relation to their 2003 class areas, buildings increased by 2.3% and transport/construction areas increased by 12% over the study time period. Coniferous forest and grasslands/open fields were the vegetated land-cover types most impacted by this urban growth, with a loss greater than 100 ha each. Expansion of the transport network negatively affected the core areas of the green infrastructure, while construction of buildings had a greater impact on habitat for species of conservation concern. Dispersal zones were affected by both types of urban growth and had the largest area converted (28 ha). Overall, there was a 9% increase in urban areas within the green infrastructure of Stockholm City. Investigation of change in habitat connectivity and availability revealed a more significant impact on probable network connectivity than habitat area for the crested tit, while by contrast common toad habitat area was more negatively affected than connectivity of the network. The impacts were due to loss, shrinkage and fragmentation of habitat patches, mainly as a result of the expansion of the transport network near the northwest and southeast outskirts of the city. Examination of the link structure of the network highlighted vulnerable or isolated parts possibly in need of strengthening or restorative efforts. Calculation of change in urban and green areas within administrative boundaries revealed that Stockholm City still has a significant proportion of green areas in comparison to urban areas, roughly 62%/38%, even if these green areas decreased slightly (-2%) as a result of urban growth over the 15-year study period.
The land-cover changes reported here are net increases in urban areas and net decreases of green areas. Yet, it is likely that there were smaller changes in the opposite direction during the study period. Indeed, during the classification process, examples of "regreening" in Stockholm City became evident. For example, an area that was previously a construction/industrial site became grassland in west Hässelby-Vällingby, buildings were removed and replaced by grassland in Nybygget in Skarpnäck and patches of barren land can be found north of Solvalla in Bromma where there were formerly buildings. There was also an example of the restoration of an edge of the green infrastructure core area in Skarpnäck as can be seen in Figure 13. Spatiotemporal land-cover analysis is effective in quantifying overall change, and classified data can be manipulated and analyzed to bring out both kinds of change. However, if an investigation requires monitoring all types of land-cover conversions, change detection techniques may be a more straightforward or practical means of highlighting changes in both directions or may be incorporated with classification techniques to provide more detailed change information, e.g., [87]. these green areas decreased slightly (-2%) as a result of urban growth over the 15-year study period.
The land-cover changes reported here are net increases in urban areas and net decreases of green areas. Yet, it is likely that there were smaller changes in the opposite direction during the study period. Indeed, during the classification process, examples of "regreening" in Stockholm City became evident. For example, an area that was previously a construction/industrial site became grassland in west Hässelby-Vällingby, buildings were removed and replaced by grassland in Nybygget in Skarpnäck and patches of barren land can be found north of Solvalla in Bromma where there were formerly buildings. There was also an example of the restoration of an edge of the green infrastructure core area in Skarpnäck as can be seen in Figure 13. Spatiotemporal land-cover analysis is effective in quantifying overall change, and classified data can be manipulated and analyzed to bring out both kinds of change. However, if an investigation requires monitoring all types of land-cover conversions, change detection techniques may be a more straightforward or practical means of highlighting changes in both directions or may be incorporated with classification techniques to provide more detailed change information, e.g., [87]. The habitat network analysis is limited to and by the city boundary, which is administrative rather than ecological in nature. Impacts on connectivity have mainly been detected on the outskirts of the study area, and for this reason, it would be interesting to shift the extent of examination to see if these impacts have more far-reaching consequences for habitat connectivity than it is possible to observe here. It should also be noted that the input parameters to the Conefor model come with some uncertainties, which will need to be further explored. In addition, the crested tit and common The habitat network analysis is limited to and by the city boundary, which is administrative rather than ecological in nature. Impacts on connectivity have mainly been detected on the outskirts of the study area, and for this reason, it would be interesting to shift the extent of examination to see if these impacts have more far-reaching consequences for habitat connectivity than it is possible to observe here. It should also be noted that the input parameters to the Conefor model come with some uncertainties, which will need to be further explored. In addition, the crested tit and common tit were selected as model species and a suite of representative and sensitive species will be needed to fulfill the target of representing the green infrastructure of the City. However, in this study, the purpose is to demonstrate the model development and novel linkage of remotely sensed land-cover change with habitat network analysis.
The monitoring of urban green infrastructure is an important part of sustainable urban development, supported by several policy documents, e.g., [88,89]. However, a recent study by Wellmann et al. [90] showed limited integration between ecology, remote sensing and urban planning, and that actionable knowledge and methods for planners and policymakers are rarely developed. This calls for method development of both planning tools and monitoring schemes. In fact, only a few studies have used medium-to high-resolution satellite data for urban green space mapping and change monitoring, e.g., [17,76]. The urban green space mapped, however, only included larger urban green spaces such as forest and parks while gardens and small green spaces were classified into other classes. For example, gardens were included in low-density built-up areas. Using Sentinel-2 data, Kopecká et al. [91] attempted to classify urban green space into 15 different classes including gardens and allotments. However, the reliability of their results is uncertain due to the lack of accuracy assessment. The highly heterogeneous urban landscape with a complex mix of impervious surfaces and vegetation requires mapping at a finer scale. As a result, very high-resolution satellite data have been exploited to classify urban green spaces in Europe and around the world, e.g., [28,29,92,93]. These studies, however, mainly focused on methods to map urban green space. No studies were found to utilize their urban green space result to further develop planning tools and monitoring schemes. This research is the first attempt to integrate high-resolution data and ecological tools for urban green space monitoring and environmental impact analysis to support planning actions. The spatiotemporal analysis of the green infrastructure and habitat network analysis could play a useful role in the local planning authority's environmental assessment and monitoring.
Stockholm City as well as Stockholm Region have developed a formal policy for biodiversity management and have been using several planning tools including green infrastructure. Consequently, a viable green infrastructure is one of the goals for urban development according to the comprehensive plan of Stockholm City [9]. Stockholm City's proposed biodiversity and related ecosystem service monitoring program [89] relies in part on locally conducted species inventories and also has planned measures such as habitat network analysis and spatiotemporal assessment of change in mapped biotopes. It mentions possible specific response measures such as installation and restoration of amphibian ponds. Results of the common toad network analysis revealed modifications to the network and provided information on where habitat has been lost and where reinstallation might be most beneficial to the network. In this way, status and change information from habitat network analysis could inform and help target restorative or conservation measures. However, there is a need for further development of these tools including the underlying habitat network analyses, and monitoring frameworks and practices. The analyses performed here are examples of practical development and implementation of habitat network and green infrastructure monitoring. The habitat networks need to be expanded to include main biodiversity components and ecosystem services, and be integrated with smart monitoring schemes to follow the changes. Lastly, the state and change of prioritized habitat networks need to be integrated into urban planning practice in order to reach sustainable urban development.

Conclusions
In this research, object-oriented SVM classifications of WorldView-2 and QuickBird-2 satellite data were used to investigate urban expansion in Stockholm City between 2003 and 2018 and the resulting environmental impacts on green areas and green infrastructure. The classifications had similar high accuracies (kappa: 0.90 and 0.92). Quantitative analysis based on the classifications has estimated broad changes, i.e., land-cover change at city landscape, down to detailed impacts, i.e., changes in terms of habitat connectivity and availability of individual habitat patches, by examining different analysis levels in Stockholm City. These are in relation to urban growth within city districts, green infrastructure and by evaluating changes in habitat connectivity and availability for selected sensitive species in the Stockholm City context. This information provides insight into where and relatively how much change occurred in different areas within the city and can inform decisions about where future urban growth should or can take place.
Between 2003 and 2018, urban areas increased and nonurban land-cover decreased by 1.7% of the city landscape area or 327 ha. Green areas decreased by 2% while urban areas increased by 4% in comparison with their 2003 aerial extent. An increase in the transport network, paved surfaces and construction sites made up 60% of the urban growth that occurred in Stockholm City between 2003 and 2018, while the remaining 40% was a net increase in buildings. This growth affected the city's green infrastructure: core areas and dispersal zones were mainly impacted by the expansion of transport and construction sites while habitat for species of conservation concern was more affected by the addition of new buildings. Expansion of the transport network and construction sites also negatively influenced connectivity in the crested tit habitat network through patch shrinkage and fragmentation mainly on the outskirts of the city. The common toad habitat network was more affected by an outright loss of patches due to the expansion of both residential areas and transport infrastructure. A continually updated spatiotemporal analysis of the green infrastructure and habitat network analysis could play a useful role in future environmental assessment by planning authorities and increase the effectiveness of restorative efforts and resources for biodiversity conservation in cities.

Common Toad
Dispersal distance Up to 2 km over optimal land-cover but only a few meters over unsuitable land-cover (Table A3) [ [96][97][98][99] Assumed maximum dispersal/movement distance in urban environment: 2 km Table A2. Resource profile for the common toad according to Mörtberg et al. [66].

Habitat Quality Biotope/Land-Cover
Reproductive habitat Optimal Wetlands, water with vegetation, streams, water-filled ditches Summer habitat Optimal Deciduous forest, mixed forest, moist grassland Suboptimal Low-density built-up areas with trees/bushes, older coniferous forest, managed grassland