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))
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