Plot filtered patch sizes in R terra

rterra

I would like to plot patch sizes above a threshold in R terra.

f=system.file("ex/elev.tif", package="terra")
r=rast(f)
x=app(r,function(x) ifelse(x< 400,NA,x))
y=patches(x)
z=tapply(cellSize(y,unit="ha")[], y[], sum)
filtered=z[z>100]

The "filtered" variable contains patches larger than 100 ha. How do I plot them in terra?

Best Answer

The filtered variable is a table of patch numbers and areas. The names of the table are the patch numbers.

To plot a map of TRUE/FALSE for those patches, do a boolean (logical) operator to see which values of y are in that set:

plot(y %in% as.numeric(names(filtered)))

or if you want a raster with those patches and their patch values (first I've made a copy so I don't stomp on y), set all that arent in the set to NA:

yf = y
yf[!yf %in% as.numeric(names(filtered))] = NA
plot(yf)

Related Solutions

Get patch size in R terra

It is best to avoid using x[] with SpatRasters as that forces all data into memory. That works fine on these examples, but may fail with your real data. In this case, you can use zonal in stead.

You example data:

library(terra)
f = system.file("ex/elev.tif", package="terra")
r = rast(f)
x = classify(r, cbind(-Inf, 400, NA))
y = patches(x)

To count the number of cells by patch

f = freq(y)
head(f)
#     layer value count
#[1,]     1     1   866
#[2,]     1     2     3
#[3,]     1     3     1
#[4,]     1     4     7
#[5,]     1     5     1
#[6,]     1     6     1

If you want the area covered by these cells, first use cellSize

z = cellSize(y,unit="km") |> zonal(y, sum)
head(z) |> round(2)
#  patches   area
#1       1 479.33
#2       2   1.66
#3       3   0.55
#4       4   3.88
#5       5   0.55
#6       6   0.55

Cells per category in R terra

If you call help(package="terra") you can bring up the help documents and read the function summaries. The function you are looking for is terra::freq

Best Answer

Related Solutions

Get patch size in R terra

Cells per category in R terra

Related Question