extract_lsm() returning wrong values

[szylinski]

I've been finding that most (not all, somehow) values returned by extract_lsm() are very wrong, when I compare them to manually measured patches in ArcMap.

Here's the code I'm running:
siteSample <- extract_lsm(MLR.r, y = sites, extract_id = sites$SiteID, what = "lsm_p_area", directions = 4)

The MLR.r raster contains 9 land cover classes of my study area, with a 10m resolution, and sites is an sp object (I have also tried an sf version of it with the same result).

This returns

> head(siteSample)
# A tibble: 6 x 7
  layer level class    id metric value extract_id
  <int> <chr> <int> <int> <chr>  <dbl>      <int>
1     1 patch     1   257 area    1.2           1
2     1 patch     1    75 area    3.68          4
3     1 patch     1   336 area    2.29          5
4     1 patch     1   336 area    2.29          8
5     1 patch     1    75 area    3.68         15
6     1 patch     1   313 area    8.31         16

Here, I know from ArcMap that the area of the site 1 patch is more like 160 hectares. Also note that sites 5 and 8 are in the same patch - they return the same, wrong value (it should be around 70).

2 or 3 sites also return the incorrect class in which they lie, so perhaps it is some projection issue, although the vast majority of sites have the correct class.

One final weird problem I ran into when troubleshooting: when I did a test run of the extract_lsm() function on a smaller cropped area of my raster and sites, some sites that returned the wrong values initially were correct, but many were still wrong.

[Nowosad]

@szylinski Could you share some version of your data with us?

[szylinski]

@Nowosad of course: landscapemetrics_test.zip. These files (a raster and a shapefile) are a small clip of my MLR.r landscape. And here's my code and output from those files:

> extract_lsm(MLRclip, y = sitesclip, extract_id = sitesclip$SiteID, what = "lsm_p_area", directions = 4)
# A tibble: 14 x 7
   layer level class    id metric  value extract_id
   <int> <chr> <int> <int> <chr>   <dbl>      <int>
 1     1 patch     1    38 area     2.93        196
 2     1 patch     1    38 area     2.93        197
 3     1 patch     1    38 area     2.93        198
 4     1 patch     1    38 area     2.93        199
 5     1 patch     1    38 area     2.93        200
 6     1 patch     1    38 area     2.93        201
 7     1 patch     1    38 area     2.93        202
 8     1 patch     1    38 area     2.93        203
 9     1 patch     1    38 area     2.93        204
10     1 patch     1    37 area   497.          205
11     1 patch     1    29 area    13.2         206
12     1 patch     1    28 area    43.5         207
13     1 patch     4   131 area     4.88        208
14     1 patch     1    12 area    20.7         209

It has correctly identified that sites 196:204 are in the same patch, but their value is wrong (should be about 445 hectares).

I also notice that site 208 is given class 4, but using raster::extract() returns the correct class 3.

[szylinski]

I've just rebuilt my code from scratch to tidy it up, and now it all seems to be working correctly... so not sure what was going wrong for me! I'm guessing this means odds are you won't be able to reproduce my problem

[mhesselbarth]

There actually was a small bug - the directions argument was not passed on correctly to the calculations of the metrics and therefore patch delineation was wrong (or still set to the default 8 at some places). Should be fixed 6df4fc3.

Please install the GitHub version (devtools::install_github("r-spatialecology/landscapemetrics) and try to rerun the code. Will also be part of the next CRAN update.

Thanks for pointing out.

library(landscapemetrics)
library(raster)
#> Loading required package: sp

MLR.r <- raster::raster("C:/Users/Maximilian/Downloads/landscapemetrics_test/MLRclip.tif")
sites <- raster::shapefile("C:/Users/Maximilian/Downloads/landscapemetrics_test/sitesclip.shp")

extract_lsm(landscape = MLR.r, y = sites, extract_id = sites$SiteID, 
            what = "lsm_p_area", directions = 4)
#> Warning: Only using 'what' argument.
#> # A tibble: 14 x 7
#>    layer level class    id metric  value extract_id
#>    <int> <chr> <int> <int> <chr>   <dbl>      <int>
#>  1     1 patch     1    38 area   497.          196
#>  2     1 patch     1    38 area   497.          197
#>  3     1 patch     1    38 area   497.          198
#>  4     1 patch     1    38 area   497.          199
#>  5     1 patch     1    38 area   497.          200
#>  6     1 patch     1    38 area   497.          201
#>  7     1 patch     1    38 area   497.          202
#>  8     1 patch     1    38 area   497.          203
#>  9     1 patch     1    38 area   497.          204
#> 10     1 patch     1    37 area    19.4         205
#> 11     1 patch     1    29 area    43.5         206
#> 12     1 patch     1    28 area    89.7         207
#> 13     1 patch     3   131 area     3.94        208
#> 14     1 patch     1    12 area    20.7         209

^{Created on 2019-09-10 by the reprex package (v0.3.0)}

[mhesselbarth]

Please re-open if necessaire