0

I have a set of points placed within one patch. I would like to measure the nearest distance between the point and the closest cell with NA values (eg. outside of my patch). Patch is defined as a cells of the same values.

Here is small scheme:

enter image description here

I am using terra::distance do measure a distance to the nearest NA value, but it seems not to work properly, and calculating the distance only for the NA cell (not from all of the other cells as I wish).

Maybe there is another easy way around, using raster or landscapemetrics approach? I would go for the simplest solution (maybe convert a point to raster first?) as I have around 50.000 points.

I have found this example but maybe there is a simpler workaround?

My working example:

library(terra)
library(raster)

# Create a small example raster
mat <- matrix(c(1, 2, 2, 2, 2,
                2, 2, 2, 2, 2,
                2, 2, 2, 2, 2,
                2, 2, 2, 2,2,
                2, 2, 2, 2, 2), byrow = T, nrow=5, ncol=5)
example_raster <- rast(nrows=5, ncols=5, vals=mat, crs = 'EPSG:3035')

# Reclassify the raster - change values '2' to NA
reclass_matrix <- matrix(c(1, NA), ncol=2, byrow=TRUE)
reclassified_raster <- classify(example_raster, reclass_matrix)

# Create a data frame with the coordinates of the point
point_df <- data.frame(x = 150, y =-80)

# Create a SpatVector from the data frame
point_vector <- vect(point_df, geom = c("x", "y"), crs = crs(reclassified_raster))

# Calculate the distance to the nearest NA,  does not work!!!
distance_to_na <- distance(reclassified_raster)     

# Plotting to visualize the results
plot(example_raster, main="example_raster")
plot(reclassified_raster, main="Reclassified Raster")
plot(distance_to_na, main="")
plot(point_vector, main="", add = T)

# Measure the distance from the point to the nearest NA
point_distance <- extract(distance_to_na, point_vector)
(point_distance)
plot(distance_to_na, main="Distance to Nearest NA")
plot(point_vector, main="", add = T)
   
3
  • 1
    The terra::distance does work, you just coded you solution incorrectly. If your reclassification matrix is matrix(c(2, NA), ncol=2, byrow=TRUE), thus recoding all your non-NA values to NA's then you will get the correct distance raster. The distance function derives distances from all NA's to all values in the raster. In your example the only value in "reclassified_raster" should be 1. You just have the problem inverted. From there you just extract the pixel value at the point and are done. For forest patches just invert the values where patches are NA and outside (current NA) are 1. Dec 6, 2023 at 19:51
  • @JeffreyEvans, thank you for this solution, it works! a stupid typo, I knew there is something inverted :-D
    – maycca
    Dec 6, 2023 at 21:36
  • @JeffreyEvans, would you like to post it as an answer? I can then accept it. Thank you
    – maycca
    Dec 7, 2023 at 3:01

1 Answer 1

0

Here is a solution converting raster patches into polygons, and subsequently into lines, using this post

library(terra)
set.seed(123)

# Create an empty raster and set all values to NA initially
raster <- rast(ext=c(401000, 402000, 4405000, 4406000), res = 30, crs="local")
values(raster) <- NA

# Define a patch area and set these cells to 2 (covering about 1/4 of the raster)
patch_cells <- cells(raster, ext(c(401000, 401500, 4405000, 4405500)))
raster[patch_cells] <- 2

# Create points that fall within the patch
points <- data.frame(
  ID = 1:4,
  x = runif(4, 401000, 401500),  # X coordinates within the patch
  y = runif(4, 4405000, 4405500) # Y coordinates within the patch
)

Test with a single point to test the results visually:

points <- data.frame(
  ID = 1,
  x = 401200,  # X coordinates within the patch
  y = 4405100 # Y coordinates within the patch
)

# convert dataframe to vector
pts <- vect(points, geom=c("x", "y"), crs=crs(raster))

# Plotting to visualize the raster and points
plot(raster, main="Raster with Patch")
plot(lns)
plot(pts, add=TRUE, col="red", pch=16)


# convert raster patches to polygons and to lines to measure the distances
lns <- as.lines(as.polygons(raster, dissolve=TRUE))

Get nearest distances:

nearest(pts, lns) |> values()

Leading to output:

 from_id from_x  from_y to_id distance
1       1 401200 4405100     1      100

enter image description here

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.