`util_extract_multibuffer` as a separate function about landscapetools

[Nowosad]

I added a PR based on your branch (it is easier to see your changes) - #36.
@bniebuhr can you create a reprex with an example or two of your new code?

[Nowosad]

I do not know the best practices. For me it is easier to understand changes when there is: a) a code comparison (PR) and b) a set of examples (reprex).

[bniebuhr]

Hi,

Sorry about the delay in woking on the examples. I made another PR with further examples and a correction of a bug - please see it here: #39
I have also included the update documentation for these two functions in the PR.

[Nowosad]

@bitbacchus @mhesselbarth @mattfrit

[bniebuhr]

Nice. I had added a few examples under @examples in the roxygen documentation but I just found a bug, so I'll work on that and commit it, and I come back with more comprehensive examples soon.

PS: just to understand the best practices, would it be better to make a PR to start the discussion then, apart from the issue?

[maelle]

Note that the NAMESPACE doesn't list util_extract_multibuffer() yet so examples using it fail as it's not available. Could you please run document()?

[mhesselbarth]

Can we close this then?

[bniebuhr]

I think there is still room for some improvement/discussion. Not sure if this should come in a separate issue, but I'll keep it here so far.

One thing is how the output of the function looks like. Testing it now, we have the following for the calculation of amount of land use classes in a classified landscape:

library(landscapetools)

# create single point
new_point = matrix(c(75,75), ncol = 2)

# show landscape and point of interest
show_landscape(classified_landscape, discrete = TRUE) +
ggplot2::geom_point(data = data.frame(x = new_point[,1], y = new_point[,2]),
                    ggplot2::aes(x = x, y = y),
                    col = "grey", size = 3)
#> Loading required package: raster
#> Loading required package: sp

# extract frequency of each pixel value within each buffer from 10 to 50 m width
util_extract_multibuffer(classified_landscape, new_point, 10, 50)
#> # A tibble: 14 x 4
#>    id          layer  freq buffer
#>    <chr>       <fct> <int>  <dbl>
#>  1 Point ID: 1 1        80     10
#>  2 Point ID: 1 2       236     10
#>  3 Point ID: 1 1       292     20
#>  4 Point ID: 1 2       964     20
#>  5 Point ID: 1 3         8     20
#>  6 Point ID: 1 1       484     30
#>  7 Point ID: 1 2      1898     30
#>  8 Point ID: 1 3       446     30
#>  9 Point ID: 1 1       825     40
#> 10 Point ID: 1 2      2810     40
#> 11 Point ID: 1 3      1389     40
#> 12 Point ID: 1 1      1362     50
#> 13 Point ID: 1 2      3750     50
#> 14 Point ID: 1 3      2748     50

^{Created on 2022-03-29 by the reprex package (v2.0.1)}

Regarding this output, I have two comments.

• First, when there is no pixels of a given class within a given buffer (e.g. no pixels of land use class 3 in the buffer with 10 m radius around this point), we do not have a new line in the output table with layer = 3 and freq = 0. Do you think having such an output would make it easier for some specific analyses (i.e. including all classes that are present in any of the buffers in the output, and filling it with 0 if it is absent from the buffer)?
• Second, I think "layer" is not really the best name for this output column (even though it is the default of raster::extract() used within the function). Maybe "class" would be better?

A second thing is how the function deals when functions are used within the buffer (instead of just counting absolute or relative frequency of pixels in each class). For example, let's say I want to take the most common land use class within each buffer size (mode function taken from here). I can do the following:

getmode <- function(v) {
    uniqv <- unique(v)
    uniqv[which.max(tabulate(match(v, uniqv)))]
}

# use a given function - e.g. mode in each buffer width
util_extract_multibuffer(classified_landscape, new_point, 10, 50, fun = getmode)
#> # A tibble: 5 x 4
#>   id          layer `function (v) \n{\n    uniqv <- unique(v)\n    uniq~` buffer
#>   <chr>       <fct>                                                 <int>  <dbl>
#> 1 Point ID: 1 2                                                         1     10
#> 2 Point ID: 1 2                                                         1     20
#> 3 Point ID: 1 2                                                         1     30
#> 4 Point ID: 1 2                                                         1     40
#> 5 Point ID: 1 2                                                         1     50

^{Created on 2022-03-29 by the reprex package (v2.0.1)}

In this case, here comes a few comments:

• Obviously there is an issue with the name of the column, we should correct that. The idea was to take the name of the function but it seems I did it wrong. Maybe we can just change to correctly take the name of the function? Or having a parameter for the user to determine this name?
• The other point is that here the column layer is the one which really represents what we want, and the other one does not mean anything. So we need to change the output when we use a function here. My suggestion is just to have Point ID, function name (e.g. getmode), value (the real value we are interested in), and buffer. Does it make sense? That would be much more appropriate for contiuous landscape maps as well (e.g. calculating some function such as the mean or weighted mean within a buffer).

I can implement all that and make a new PR but any thoughts are welcome.

[bniebuhr]

Just to complement, I guess when I deal with the last point, it will be more meaningful to use that for continuous maps as well, see example below:

library(landscapetools)

# create single point
new_point = matrix(c(75,75), ncol = 2)

# show landscape and point of interest
show_landscape(fractal_landscape) +
ggplot2::geom_point(data = data.frame(x = new_point[,1], y = new_point[,2]),
                    ggplot2::aes(x = x, y = y),
                    col = "grey", size = 3)
#> Loading required package: raster
#> Loading required package: sp

# use a given function - e.g. median in each buffer width
util_extract_multibuffer(fractal_landscape, new_point, 10, 50, fun = "mean")
#> # A tibble: 5 x 4
#>   id          layer              mean buffer
#>   <chr>       <fct>             <int>  <dbl>
#> 1 Point ID: 1 0.466132288392125     1     10
#> 2 Point ID: 1 0.475467530484216     1     20
#> 3 Point ID: 1 0.520185362184868     1     30
#> 4 Point ID: 1 0.547060865215994     1     40
#> 5 Point ID: 1 0.561445545416456     1     50

^{Created on 2022-03-29 by the reprex package (v2.0.1)}

What we want is the value (currently in the layer column).

[Nowosad]

Hi @bniebuhr, my thinking about your ideas:

1. I think it would be nice to have 0 for non-existing classes (for the given buffer)
2. The class name is fine (for that example)
3. Regarding the use of custom functions -- it is possible to use an anonymous function here. Therefore, how about just naming this column with generic fun?
4. Custom functions may have different outputs compared to the default fun = NULL... Could you make a list of the most common (in your opinion) use of the custom functions? Then we could try to think about the best solution...