Code presentation of the landscapemetrics R package within the Seminar of the EAGLE module ‘Spatial Modeling and Prediction’

Theory of Landscape metrics

Ecology and landscape ecology mainly studies interactions between organisms and their environment. In a wider sense, the environment perceived by the organism can be defined as a landscape. This can spread from huge areas for, e.g., mobile, large-bodied mammals to much smaller areas for, e.g., immobile insects and always depends on the research question.

Landscape metrics are tools to characterize a landscape. This includes mainly describing the composition and configuration of a landscape. While the composition basically describes how much of the landscape is occupied by a certain land cover type, the configuration mainly describes the spatial arrangement of the land cover types. The basic idea of landscape metrics is to condense as much information as possible into single number.

Landscapes are dynamic systems and Humans affects them continuously. Depending on intensive human effects, pressure is increased on landscapes. Consequently, landscapes are altered over time and there are negative effects of pressures on landscape and species living in the area. The negative effects are especially vulnerable more intense to the human effected landscapes. In these landscapes, fragmentation increased. Habitats have been damaged. Depending on these effects material flow and transactions of the species are limited. Landscape ecology investigates landscape structure and changes in the landscape. Change expresses any modification occurring in the landscape over time. Landscape structure evaluates land mosaic as measure, number, size and shape. (Gökyer, 2013)

landscapemetrics Package

Levels of landscape metrics

Landscape metrics can be calculated for three different levels (or “scales”). Each level contains information about different aspects of the landscape. Level of interest depends largely on the research question, and often a combination of several levels can be useful. Similar, also the decision of which metric to use depends on the research question asked.

A patch is defined as neighboring cells belonging to the same class i. Hereby, landscapemetrics uses the 8-neighbors rule (Queen’s case) to identify patches. Patch level metrics are calculated for every patch in the landscape, regardless of the class that the patch belongs to. The output will match the number of patches present in the landscape. These metrics are also often the basis for metrics of the other two levels.

Class level metrics summaries all patches belonging to one class i. These metrics can be either the “distribution” of patch level metrics of all patches of class i (e.g., the mean) or consider only patches of class i for the calculations of the metric. Regardless of the mathematical background, the output will always match the number of classes present. Class level metrics are suitable to describe the composition and configuration of the landscape.

Landscape level metrics summaries the whole landscape into one value. This can either be done by summarising metrics of lower levels or calculating a metric including all patches and classes. Following, the output will always be just one number. Landscape level metrics are suitable for condensing information about the landscape into just one value.

Example

Changes in the Composition and Configuration of a landscape at species sampling locations. Example using National Landcover Database (NLCD) of the US and species locations data of the endangered Ruffed Grouse in King County, Washington State from gbif.

Code

load or install necessary libraries

#load libraries
library(rnaturalearth)
library(sf)
library(FedData)
library(raster)
library(exactextractr)
library(dplyr)
library(terra)
library(dismo)
library(landscapemetrics)
library(rasterVis)
library(gridExtra)
library(animation)
library(ggplot2)

Area of Intrest

get area of intrest, in this case King County, in Washington State

#calling our state and county of interest 
state = c("Washington")
county <- c('King')

#getting our county and state
Wa <- getData("GADM",country="USA",level=2)
WA <- Wa[Wa$NAME_1 %in% state,]

#get the desired county and reproject it to nlcd crs
Kg<- st_as_sf(WA[WA$NAME_2 %in% county,])%>%
    st_transform(crs = "+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=23 +lon_0=-96 
    +x_0=0 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs") 

#reproject to WGS84 for species download in gbif
Kg_p <- st_transform(Kg, crs = "+proj=longlat +datum=WGS84")

NLCD-data

download NLCD-data

#get nlcd data (national land cover database of the US)
#download nlcd data via R, apparently just works or seems for gdam Level2 data,
#e.g states and counties 

#for 2019, more years are available 
nlcd01 <- get_nlcd(
  Kg,
  year = 2001,
  label = "WA",
  force.redo = T
)

nlcd19 <- get_nlcd(
  Kg,
  year = 2019,
  label = "WA",
  force.redo = T
)

#ratify data, that it's recognized to be categorical data for the later use in
#landscapemetrics
nlcd01 <- ratify(nlcd01)
nlcd19 <- ratify(nlcd19)
s <- stack(nlcd01, nlcd19)

Visualizung

visualize NLCD-data

#visualizing nlcd data with legend
#get legend from nlcd data set
legend<-pal_nlcd()
legend
## # A tibble: 20 × 4
##       ID Class                        Color   Description                       
##    <dbl> <chr>                        <chr>   <chr>                             
##  1    11 Open Water                   #5475A8 Areas of open water, generally wi…
##  2    12 Perennial Ice/Snow           #FFFFFF Areas characterized by a perennia…
##  3    21 Developed, Open Space        #E8D1D1 Areas with a mixture of some cons…
##  4    22 Developed, Low Intensity     #E29E8C Areas with a mixture of construct…
##  5    23 Developed, Medium Intensity  #ff0000 Areas with a mixture of construct…
##  6    24 Developed High Intensity     #B50000 Highly developed areas where peop…
##  7    31 Barren Land (Rock/Sand/Clay) #D2CDC0 Areas of bedrock, desert pavement…
##  8    41 Deciduous Forest             #85C77E Areas dominated by trees generall…
##  9    42 Evergreen Forest             #38814E Areas dominated by trees generall…
## 10    43 Mixed Forest                 #D4E7B0 Areas dominated by trees generall…
## 11    51 Dwarf Scrub                  #AF963C Alaska only areas dominated by sh…
## 12    52 Shrub/Scrub                  #DCCA8F Areas dominated by shrubs; less t…
## 13    71 Grassland/Herbaceous         #FDE9AA Areas dominated by gramanoid or h…
## 14    72 Sedge/Herbaceous             #D1D182 Alaska only areas dominated by se…
## 15    73 Lichens                      #A3CC51 Alaska only areas dominated by fr…
## 16    74 Moss                         #82BA9E Alaska only areas dominated by mo…
## 17    81 Pasture/Hay                  #FBF65D Areas of grasses, legumes, or gra…
## 18    82 Cultivated Crops             #CA9146 Areas used for the production of …
## 19    90 Woody Wetlands               #C8E6F8 Areas where forest or shrubland v…
## 20    95 Emergent Herbaceous Wetlands #64B3D5 Areas where perennial herbaceous …
#get values of the raster classes 
vals<-unique(nlcd01)

#and apply corresbonding names in legend to them
df<-legend[legend$ID %in% vals,]


# create a custom legend:
myKey <- list(rectangles=list(col = df$Color),
              text=list(lab=df$Class),
              space='left',
              columns=1,
              size=2,
              cex=.6)


# levelplot data 
levelplot(nlcd01, att='ID', 
          col.regions=df$Color,
          par.settings = list(axis.line = list(col = "transparent"), 
                              strip.background = list(col = 'transparent'), 
                              strip.border = list(col = 'transparent')), 
          scales = list(col = "transparent"),
          colorkey=F,
          key=myKey,
          main="Landcover 2001")

levelplot(nlcd19, att='ID', 
                col.regions=df$Color,
                par.settings = list(axis.line = list(col = "transparent"), 
                                    strip.background = list(col = 'transparent'), 
                                    strip.border = list(col = 'transparent')), 
                scales = list(col = "transparent"),
                colorkey=F,
                key=myKey,
                 main="Landcover 2019")

GIF

animation of land cover change between 2001 and 2019

# Make a list of years in string format for titles
years.list<-list("2001",  "2019")

#crate animation of land cover change during the years 
saveGIF(
  {
    for(i in c(1:nlayers(s))){
      rat<-ratify(s[[i]])
      a<-levelplot(rat, att='ID', 
                   col.regions=df$Color,
                   par.settings = list(axis.line = list(col = "transparent"), 
                                       strip.background = list(col = 'transparent'), 
                                       strip.border = list(col = 'transparent')), 
                   scales = list(col = "transparent"),
                   main=paste0("King County land cover ", years.list[[i]]),
                   colorkey=F,
      )
      print(a)
    }
  }, interval=0.8, movie.name="szoom.gif", ani.width = 600)
## [1] TRUE

GBIF-species data

get locations of spotted Ruffed Grouse in King County

#get species data from gbif
#download species data from gbif, in that case endangered Ruffed Grouse 
#in King county, washington state
spec_occ <- gbif("Bonasa", "umbellus", 
                 ext = extent(Kg_p), 
                 sp = TRUE, removeZeros = TRUE,
                 start = 0,
                 end = 500)


#assign a projection to the species data
proj4string(spec_occ) <- CRS("+proj=longlat +datum=WGS84")

#reproject to crs of nlcd data
spec <- spTransform(spec_occ, CRS("+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=23 
                         +lon_0=-96 +x_0=0 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 
                         +units=m +no_defs"))

#check for complete cases within study area
spec <- spec[complete.cases(extract(nlcd01, spec)), ]

#have look at the saml´pling locations 
plot(nlcd01, main="Landcover 2001 and Species Occurance")
plot(spec, add=T)

landscapemetrics

look at the landscape composition and configuration of the landcover data, and differences between years

#Extract raster values at sampling locations
spec@data <- data.frame(spec@data, extract(nlcd01, spec)) 
spec@data <- data.frame(spec@data, extract(nlcd19, spec)) 


#What land cover type is assigned to the most sampling units? 
table(spec@data$extract.nlcd01..spec.)
## 
## 11 21 22 23 24 41 42 43 52 71 81 90 
## 13 73 23 18 40  7 55 37  6  3  7 23
table(spec@data$extract.nlcd19..spec.)
## 
## 11 21 22 23 24 41 42 43 52 71 81 90 
## 13 92 24 18 40  7 52 35  4  3  7 10
#check if raster is suitable for use in landscapemetrics 
landscapemetrics::check_landscape(nlcd01)
##   layer       crs units   class n_classes OK
## 1     1 projected     m integer        16  ✔
landscapemetrics::check_landscape(nlcd19)
##   layer       crs units   class n_classes OK
## 1     1 projected     m integer        16  ✔
#have a look at the available metrics for landscape level, possible to filter for
#other levels
landscapemetrics::list_lsm(level = "landscape", type = "diversity metric")
## # A tibble: 9 × 5
##   metric name                               type             level function_name
##   <chr>  <chr>                              <chr>            <chr> <chr>        
## 1 msidi  modified simpson's diversity index diversity metric land… lsm_l_msidi  
## 2 msiei  modified simpson's evenness index  diversity metric land… lsm_l_msiei  
## 3 pr     patch richness                     diversity metric land… lsm_l_pr     
## 4 prd    patch richness density             diversity metric land… lsm_l_prd    
## 5 rpr    relative patch richness            diversity metric land… lsm_l_rpr    
## 6 shdi   shannon's diversity index          diversity metric land… lsm_l_shdi   
## 7 shei   shannon's evenness index           diversity metric land… lsm_l_shei   
## 8 sidi   simpson's diversity index          diversity metric land… lsm_l_sidi   
## 9 siei   simspon's evenness index           diversity metric land… lsm_l_siei
#check available metrics filtered by metric
landscapemetrics::list_lsm(metric = "area")
## # A tibble: 7 × 5
##   metric  name       type                 level     function_name
##   <chr>   <chr>      <chr>                <chr>     <chr>        
## 1 area    patch area area and edge metric patch     lsm_p_area   
## 2 area_cv patch area area and edge metric class     lsm_c_area_cv
## 3 area_mn patch area area and edge metric class     lsm_c_area_mn
## 4 area_sd patch area area and edge metric class     lsm_c_area_sd
## 5 area_cv patch area area and edge metric landscape lsm_l_area_cv
## 6 area_mn patch area area and edge metric landscape lsm_l_area_mn
## 7 area_sd patch area area and edge metric landscape lsm_l_area_sd
#filter by level and metric
landscapemetrics::list_lsm(level = c("class", "landscape"), type = "aggregation metric")
## # A tibble: 29 × 5
##    metric   name                                  type       level function_name
##    <chr>    <chr>                                 <chr>      <chr> <chr>        
##  1 ai       aggregation index                     aggregati… class lsm_c_ai     
##  2 clumpy   clumpiness index                      aggregati… class lsm_c_clumpy 
##  3 cohesion patch cohesion index                  aggregati… class lsm_c_cohesi…
##  4 division division index                        aggregati… class lsm_c_divisi…
##  5 enn_cv   euclidean nearest neighbor distance   aggregati… class lsm_c_enn_cv 
##  6 enn_mn   euclidean nearest neighbor distance   aggregati… class lsm_c_enn_mn 
##  7 enn_sd   euclidean nearest neighbor distance   aggregati… class lsm_c_enn_sd 
##  8 iji      interspersion and juxtaposition index aggregati… class lsm_c_iji    
##  9 lsi      landscape shape index                 aggregati… class lsm_c_lsi    
## 10 mesh     effective mesh size                   aggregati… class lsm_c_mesh   
## # ℹ 19 more rows
#get percentage of class occupied in the landscape in 2001 and 2019, bind them
percentage_class01 <- lsm_c_pland(landscape = nlcd01)
percentage_class19 <- lsm_c_pland(landscape = nlcd19)
names(percentage_class01) <- c("layer01", "level01", "class01", "id01", 
                               "metric01", "value01")
percentage_class <- bind_cols(percentage_class01, percentage_class19)

#calculate the difference in percentage of classes occupied in the landscape
percentage_class$dif <- percentage_class$value-percentage_class$value01
percentage_class
## # A tibble: 16 × 13
##    layer01 level01 class01  id01 metric01 value01 layer level class    id metric
##      <int> <chr>     <int> <int> <chr>      <dbl> <int> <chr> <int> <int> <chr> 
##  1       1 class        11    NA pland      7.10      1 class    11    NA pland 
##  2       1 class        12    NA pland      0.116     1 class    12    NA pland 
##  3       1 class        21    NA pland      5.19      1 class    21    NA pland 
##  4       1 class        22    NA pland      6.42      1 class    22    NA pland 
##  5       1 class        23    NA pland      4.76      1 class    23    NA pland 
##  6       1 class        24    NA pland      1.99      1 class    24    NA pland 
##  7       1 class        31    NA pland      1.37      1 class    31    NA pland 
##  8       1 class        41    NA pland      1.84      1 class    41    NA pland 
##  9       1 class        42    NA pland     45.7       1 class    42    NA pland 
## 10       1 class        43    NA pland      8.90      1 class    43    NA pland 
## 11       1 class        52    NA pland      8.93      1 class    52    NA pland 
## 12       1 class        71    NA pland      2.86      1 class    71    NA pland 
## 13       1 class        81    NA pland      2.60      1 class    81    NA pland 
## 14       1 class        82    NA pland      0.113     1 class    82    NA pland 
## 15       1 class        90    NA pland      1.68      1 class    90    NA pland 
## 16       1 class        95    NA pland      0.427     1 class    95    NA pland 
## # ℹ 2 more variables: value <dbl>, dif <dbl>
#Visualize patches 
# get patches for each class and year 
cc_nlcd01 <- landscapemetrics::get_patches(nlcd01, directions = 8)
cc_nlcd19 <- landscapemetrics::get_patches(nlcd19, directions = 8)

#plot them selected class to see difference
plot(cc_nlcd01$layer_1$class_42, main="Evergreen Forest Patch 2001")

plot(cc_nlcd19$layer_1$class_42, main="Evergreen Forest Patch 2019")

landscapemetrics at sampling locations

extract patch area for all classes at species locations

#extract patch area for all classes at species locations
patch_size_sp01 <- extract_lsm(landscape = nlcd01, y = spec, what = "lsm_p_area")
patch_size_sp19 <- extract_lsm(landscape = nlcd19, y = spec, what = "lsm_p_area")
patch_size <- cbind(patch_size_sp01, patch_size_sp19) 

#change names 
names(patch_size) <- c("layer01", "level01", "class01", "id01", "metric01", "value01", 
"extract_id01", "layer19", "level19", "class19", "id19", "metric19","value19", "extract_id19")

## get patch size for one class(Mixed Forest), we set all area of class != 43 to 0:
patch_size_sp_forest01 <- dplyr::mutate(patch_size_sp01, 
                                      value = dplyr::case_when(class == 43 ~ value, 
                                                               class != 43~ 0))
patch_size_sp_forest19 <- dplyr::mutate(patch_size_sp19, 
                                      value = dplyr::case_when(class == 43 ~ value, 
                                                               class != 43~ 0))

## a.d data to species data:
spec@data$ForestPatchSize01 <- patch_size_sp_forest01$value
spec@data$ForestPatchSize19 <- patch_size_sp_forest19$value


#visualize patches to compare them 
b1 <- bubble(spec, "ForestPatchSize01", fill = FALSE, 
      key.entries = as.numeric(names(table(spec@data$ForestPatchSize01))))
b2 <- bubble(spec, "ForestPatchSize19", fill = FALSE, 
       key.entries = as.numeric(names(table(spec@data$ForestPatchSize19))))
plot(b1)

plot(b2)

Calculate metrics for species sampling locations

land_met01 <- sample_lsm(nlcd01, 
                          what = c("lsm_l_contag", #Contagion
                                    "lsm_l_pd", #Patch density
                                    "lsm_l_ed", #Edge density
                                    "lsm_l_lsi", #landscape shape index
                                    "lsm_l_lpi", #Largest patch index
                                    "lsm_l_pr"),#Patch richness
                          shape = "square",
                          y = spec, 
                          size = 500,
                          directions = 8)
#This gives the layer a name associated with the landscape so we can bind them!
land_met01$layer <- "2001" 

#same for 2019
land_met19 <- sample_lsm(nlcd19, 
                          what = c("lsm_l_contag", #Contagion
                                   "lsm_l_pd", #Patch density
                                   "lsm_l_ed", #Edge density
                                   "lsm_l_lsi", #landscape shape index
                                   "lsm_l_lpi", #Largest patch index
                                   "lsm_l_pr"),#Patch richness
                          shape = "square",
                          y = spec, 
                          size = 500,
                          directions = 8)
land_met19$layer <- "2019"

#crate a list of the data tables and bind them
allmetrics1 <- list(land_met01, land_met19) 
metric_table1 <- do.call(rbind, allmetrics1) 
allmetrics1
## [[1]]
## # A tibble: 1,830 × 8
##    layer level     class    id metric  value plot_id percentage_inside
##    <chr> <chr>     <int> <int> <chr>   <dbl>   <int>             <dbl>
##  1 2001  landscape    NA    NA contag  51.8        1              98.0
##  2 2001  landscape    NA    NA ed     229.         1              98.0
##  3 2001  landscape    NA    NA lpi     19.9        1              98.0
##  4 2001  landscape    NA    NA lsi      6.67       1              98.0
##  5 2001  landscape    NA    NA pd      70.4        1              98.0
##  6 2001  landscape    NA    NA pr      14          1              98.0
##  7 2001  landscape    NA    NA contag  36.3        2              98.0
##  8 2001  landscape    NA    NA ed     281.         2              98.0
##  9 2001  landscape    NA    NA lpi     10.5        2              98.0
## 10 2001  landscape    NA    NA lsi      7.95       2              98.0
## # ℹ 1,820 more rows
## 
## [[2]]
## # A tibble: 1,830 × 8
##    layer level     class    id metric  value plot_id percentage_inside
##    <chr> <chr>     <int> <int> <chr>   <dbl>   <int>             <dbl>
##  1 2019  landscape    NA    NA contag  51.7        1              98.0
##  2 2019  landscape    NA    NA ed     229.         1              98.0
##  3 2019  landscape    NA    NA lpi     19.9        1              98.0
##  4 2019  landscape    NA    NA lsi      6.66       1              98.0
##  5 2019  landscape    NA    NA pd      67.3        1              98.0
##  6 2019  landscape    NA    NA pr      14          1              98.0
##  7 2019  landscape    NA    NA contag  36.0        2              98.0
##  8 2019  landscape    NA    NA ed     289.         2              98.0
##  9 2019  landscape    NA    NA lpi      8.08       2              98.0
## 10 2019  landscape    NA    NA lsi      8.14       2              98.0
## # ℹ 1,820 more rows
#visualize in ggplot 
p <- ggplot(data = metric_table1, mapping = aes(x = as.factor(layer), y = value))
p + geom_bar(stat = "identity") + 
  facet_wrap(~metric, scales = "free") +
  xlab("Landscape") +
  ggtitle("Landscape level metrics for 2001 and 2019 landscapes")

Implications

Fragmentation refers to the phenomena of a large patch/habitat/land cover type transforming into a number of smaller patches/habitats/land cover types that are not as spatially connected to each other, if at all (as they were in their former state).

Aggregation

  • Patch density (‘pd’) is a metric which describes fragmentation. The patch density metric value increases as the landscape gets more fragmented or patchy. As such, it appears that the 2019 landscape is more fragmented than the 2001 landscape at the sampling locations.

  • Landscape shape index (a ratio between the actual landscape edge length and the hypothetical minimum edge length. The LSI metric value increases as edge length increases, or as the patches become less compact)

Area and Edge

  • Edge density (the equal of all edges in the landscape, describes the configuration of the landscape (ie. highly aggregated classes will result in a low edge density).);

  • Largest patch index (the percentage of the landscape covered by the corresponding largest patch of each class. It is a measure of dominance)

Diversity - Patch richness (measures the number of classes in a landscape)

Landscape change over time based on class level

land_mets01 <- sample_lsm(nlcd01, 
                          what = c(
                            "lsm_c_pd", #Patch density
                            "lsm_c_ed", #Edge density
                            "lsm_c_lsi", #landscape shape index
                            "lsm_c_lpi" #Largest patch index
                          ),
                          shape = "square",
                          y = spec, 
                          size = 500,
                          directions = 8)
land_mets01$layer <- "2001"

land_mets19 <- sample_lsm(nlcd19, 
                          what = c(
                            "lsm_c_pd", #Patch density
                            "lsm_c_ed", #Edge density
                            "lsm_c_lsi", #landscape shape index
                            "lsm_c_lpi" #Largest patch index
                          ),
                          shape = "square",
                          y = spec, 
                          size = 500,
                          directions = 8)
land_mets19$layer <- "2019"

#creates list of the data tables and bind them
allmetrics2 <- list(land_mets19, land_mets01) 
metric_table2 <- do.call(rbind, allmetrics2) 
allmetrics2
## [[1]]
## # A tibble: 10,448 × 8
##    layer level class    id metric  value plot_id percentage_inside
##    <chr> <chr> <int> <int> <chr>   <dbl>   <int>             <dbl>
##  1 2019  class    11    NA ed       2.45       1              98.0
##  2 2019  class    21    NA ed      61.2        1              98.0
##  3 2019  class    22    NA ed      34.3        1              98.0
##  4 2019  class    23    NA ed      19.0        1              98.0
##  5 2019  class    24    NA ed       2.45       1              98.0
##  6 2019  class    31    NA ed       5.51       1              98.0
##  7 2019  class    41    NA ed     125.         1              98.0
##  8 2019  class    42    NA ed      23.9        1              98.0
##  9 2019  class    43    NA ed     155.         1              98.0
## 10 2019  class    52    NA ed      10.1        1              98.0
## # ℹ 10,438 more rows
## 
## [[2]]
## # A tibble: 10,152 × 8
##    layer level class    id metric  value plot_id percentage_inside
##    <chr> <chr> <int> <int> <chr>   <dbl>   <int>             <dbl>
##  1 2001  class    11    NA ed       2.45       1              98.0
##  2 2001  class    21    NA ed      61.2        1              98.0
##  3 2001  class    22    NA ed      33.7        1              98.0
##  4 2001  class    23    NA ed      18.4        1              98.0
##  5 2001  class    24    NA ed       1.22       1              98.0
##  6 2001  class    31    NA ed       5.51       1              98.0
##  7 2001  class    41    NA ed     125.         1              98.0
##  8 2001  class    42    NA ed      23.0        1              98.0
##  9 2001  class    43    NA ed     156.         1              98.0
## 10 2001  class    52    NA ed      13.5        1              98.0
## # ℹ 10,142 more rows
p <- ggplot(data = metric_table2, mapping = aes(fill =layer, x = factor(class), y = value))
p + geom_bar(position = "dodge", stat = "identity") +
  facet_wrap(~metric, scales = "free") +
  xlab("Class") +
  scale_fill_grey(start=0.8, end=0.2) +
  ggtitle("Landscape change at sampling locations: Between 2001 and 2019")

Percentage change of metrics

#calculate percentage 
metric_table_wide <- tidyr::spread(metric_table2, layer, value)
metric_table_wide <- metric_table_wide %>%
  mutate(percent_change = (`2019`-`2001`)/`2001`) 

#summarised percentage change 
grouped <- metric_table_wide[!is.na(metric_table_wide$percent_change),] %>%
  group_by(class, metric) %>%
  summarise(percentage_change = mean(percent_change))
grouped
## # A tibble: 60 × 3
## # Groups:   class [15]
##    class metric percentage_change
##    <int> <chr>              <dbl>
##  1    11 ed                0.0157
##  2    11 lpi              -0.0343
##  3    11 lsi               0.0396
##  4    11 pd                0.172 
##  5    21 ed                0.814 
##  6    21 lpi               1.61  
##  7    21 lsi               0.0658
##  8    21 pd                0.0433
##  9    22 ed                0.0452
## 10    22 lpi              -0.0304
## # ℹ 50 more rows

Implications

following generalized inferences can be made:

  • Patch density increased for more classes, then it decreased

  • Edge density decreased for most classes or didn’t changed

  • Landscape shape index slightly increased

  • Largest patch index slightly decreased or didn’t changed

The combination of these metrics infers that between 2001 and 2019, the landscape at the Ruffed Grouse locations became a bit more fragmented, the landscape became little less aggregated by class, the patches became less compact and the landscape also became little less dominated by a particular land cover type in 2019.

References

  • Hesselbarth, M.H.K., Sciaini, M., With, K.A., Wiegand, K., Nowosad, J. 2019. landscapemetrics: an open‐source R tool to calculate landscape metrics. Ecography, 42: 1648-1657 (v0.0)

  • Hesselbarth, M, Wagner Helena. Worked Example. https://bookdown.org/hhwagner1/LandGenCourse_book/WE_2.html

  • Gökyer, Ercan: Understanding Landscape Structure Using Landscape Metrics, July 1st, 2013.

  • Nowosad J., TF Stepinski. 2019. Information theory as a consistent framework for quantification and classification of landscape patterns.

  • White, Emma.2022. Landscape Ecology – Prac 1. https://rpubs.com/whitee2/882620