How do I replicate this raster script using terra (specifically struggling with the mask function)?

Viewed 302

I want to extract some application rates, total application and the area from a raster file containing pesticide use. I have many files where I need to do this and raster is to slow so I need to do it in terra instead but struggling to replicate the mask function.

Here is the code in raster (copied from here: Using raster to calculate the mean application and total application of pesticides, but numbers not adding up):

library(tidyverse)

## With raster ----
data(wrld_simpl)
r <- raster::raster("https://raw.github.com/hansronald/Pesticide-data/master/APR_Soybean_Glyphosate_2015_L.tif")
r <- raster::clamp(r, lower=0, useValues=FALSE)

# area in ha
a <- raster::area(r) * 100

## Get the mean application rate
mean_app <- raster::extract(r, wrld_simpl, fun = mean, na.rm = TRUE)
rtot  <- r * a

## Get the total application for each country
tot_app  <- raster::extract(rtot, wrld_simpl, fun = sum, na.rm = TRUE)

## Get the total area for each country
rarea <- mask(a, r)
tot_area <- raster::extract(rarea, wrld_simpl, fun = sum, na.rm = TRUE)

In terra

## Terra ----

data(wrld_simpl)
wrld_simpl = vect(wrld_simpl)

r <- terra::rast("https://raw.github.com/hansronald/Pesticide-data/master/APR_Soybean_Glyphosate_2015_L.tif")
r <- terra::clamp(r, lower=0, values=FALSE)

# area in ha
a <- terra::area(r) * 100

mean_app <- terra::extract(r, wrld_simpl, fun = mean, na.rm = TRUE)
rtot  <- r * a
tot_app  <- terra::extract(rtot, wrld_simpl, fun = sum, na.rm = TRUE)

rarea <- terra::mask(a, r)
tot_area <- terra::extract(rarea, wrld_simpl, fun = sum, na.rm = TRUE)

But I get this error when I try to use mask:

Error in (function (classes, fdef, mtable)  : 
  unable to find an inherited method for function ‘mask’ for signature ‘"numeric", "SpatRaster"’

Any ideas how to solve this?

1 Answers

The area method in terra is a bit different than in raster. In terra it can return either a SpatRaster (area for each cell) or a number (the total area). In your case a is the total area, and there should be no surprise that you cannot use it with mask (it is of course good to inspect the objects at each step.)

You can get what you want by using sum=FALSE. And by adding mask=TRUE you can skip the mask step altogether.

a <- terra::area(r, sum=FALSE, mask=TRUE) / 10000

Also note that the unit is m2.

Related