0

I'm trying to calculate Topographic Position Index (TPI) for 177 points of interest. I have their coordinates stored in a data.frame and elevation in a raster of 7.5 arc sec spatial resolution. And the TPIs I'm calculating is basically: the elevation of point of interest minus the average elevation of its surrounding cells, then the intermediate result is divided by the spatial resolution of the raster. (resolution(dem)) to account for differences in the spatial scale of the DEM and the TPI values.

And since studies usually calculate two TPIs (a small scale + a large scale), I'm also using two windows, where in the small one the surrounding 55 cells are used, and in the large one the surrounding 1010 cells are used.

I am getting this error message

Error in .focal_fun(v, w, as.integer(c(tr$nrows[1] + addr, nc)), runfun,  :    
   Evaluation error: could not find function "resolution"

Code:

library(raster)
library(sp)
library(terra)
library(haven)

# Make a dataframe with longitudes and latitudes
df <- data.frame(lon = coords$longitude, lat = coords$latitude)

# Convert the dataframe to a SpatialPointsDataFrame
coordinates(df) <- c("lon", "lat")
proj4string(df) <- CRS("+proj=longlat +datum=WGS84")

# Extract the elevation values at the points
elev <- extract(dem, df)

# Define the scales for TPI calculation
scales <- c(5, 10)

# Loop over the scales and calculate TPI
tpi_list <- list()
for (scale in scales) {
  # Define the size of the moving window
  win_size <- scale * 5
  
  # Calculate TPI
  tpi <- focal(dem, w = matrix(1, win_size, win_size), fun = function(x) {
    (elev - mean(x)) / resolution(dem) * 5
  })

  # Extract TPI values at the points
  tpi_vals <- extract(tpi, df)
  
  # Store the TPI values in a list
  tpi_list[[as.character(scale)]] <- tpi_vals
}

#Error in .focal_fun(v, w, as.integer(c(tr$nrows[1] + addr, nc)), runfun,  :    Evaluation error: could not find function "resolution"

# Combine the TPI values for different scales into a dataframe
tpi_df <- data.frame(tpi_list, row.names = rownames(df))
Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63
  • 1
    clearly the error occurs in `focal`. So please change your example and make it *minimal*. Remove the loop, and anything else that you do not need to produce the error. Also make it self contained and reproducible by including example data created with code or using a file that ships with R (e.g the elevation data that come with "terra", see `?rast` – Robert Hijmans Feb 04 '23 at 21:19

1 Answers1

1

The error you get is clear:

Evaluation error: could not find function "resolution"

You are using a function resolution, but R does not know about that function. It does not exist in the current workspace. I suppose you were looking for res.

Here is a working example.

library(terra)
dem <- rast(system.file("ex/elev.tif", package="terra"))
df <- data.frame(lon=c(5.9, 6.0, 6.2), lat=c(49.9, 49.6, 49.7))

tpifun <- \(x, f) x[f] - mean(x[-f], na.rm=TRUE)

scales <- c(5, 11)
tpilst <- vector("list", length(scales))
for (i in seq_along(scales)) {
    win_size <- scales[i] * 5
    mid <- ceiling(win_size^2 / 2)
    tpi <- focal(dem, w=win_size, fun=tpifun, f=mid, wopt=list(names="tpi"))
    tpilst[[i]] <- data.frame(scale=scales[i], extract(tpi, df))
}

tpi <- do.call(rbind, tpilst)
tpi$tpi <- tpi$tpi / (mean(res(dem)) * 5)

tpi
#  scale ID        tpi
#1     5  1 -1330.6184
#2     5  2   340.1538
#3     5  3  -135.3077
#4    11  1  -585.7952
#5    11  2   255.3344
#6    11  3   292.0155

A couple of things:

res returns two numbers, the x and y resolution. In the example above I take the mean.

What you were doing in the function supplied to focal is not possible. You supplied a data.frame with elevation data for a few points. How can focal understand what that is all about? Instead, you can compute the TPI for each cell and extract these values.

You cannot use a value of 10 for scale because the weights matrix must have odd size. Otherwise it is not clear how it should be centered on the focal cell.

You say that if scale is 5, the "surrounding 55 cells are used". But that is not the case. The number of surrounding cells used is 624.

scale <- 5
(scale * 5)^2 - 1
#[1] 624
Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63
  • Thank you for using a system file as an example. – idontreallyknowhehe Feb 05 '23 at 00:45
  • I got this error message: ``` Error in .local(x, ...) : is.matrix(w) is not TRUE ``` So I changed my "x" in the focal function from the DEM to a matrix: dem_matrix <- as.matrix(dem) . Then, I got a new error message: Error in (function (classes, fdef, mtable) : unable to find an inherited method for function ‘focal’ for signature ‘"matrix"’ . Do you know how I can solve it? Thank you! – idontreallyknowhehe Feb 05 '23 at 03:06
  • You can copy and paste my example and see that it works. So all you need to do is replace the filename and data.frame with your values. To not get confused, do not load other packages such as "raster" – Robert Hijmans Feb 05 '23 at 18:38
  • I see. I think it's the DEM's issue. I read the DEM using the "raster" package rather than "terra". Thank you! – idontreallyknowhehe Feb 05 '23 at 23:49
  • Well, that is a mistake... – Robert Hijmans Feb 06 '23 at 06:13