Calculate zonal statistics from rasters in a zip file

Imagine you have a zipped folder with a bunch of raster maps containing variables (e.g. the ones you can download from WorldClim or from CliMond), and you need to calculate zonal statistics from each of these rasters. The zonaFromZip function, provided below, automates this process without the need to unzip the folders. It extracts one raster at a time from the .zip, imports it to R, calculates zonal statistics for your zones raster map (the ‘mean’ function is used by default, but you can provide any other argument accepted by the zonal function of the R raster package), and then deletes the unzipped file before unzipping the next one, therefore requiring minimal disk space.

zonalFromZip <- function (zip.file, zones.rast, rast.file.ext = ".tif", aux.file.ext = NULL, delete.unzipped = TRUE, verbosity = 2, ...)
  # zip.file: path to the zip containing the raster maps to calculate zonal stats from
  # zones.rast: raster map (in your R workspace) containing the spatial units to calculate zonal stats to
  # rast.file.ext: file extension for the raster maps on disk
  # aux.file.ext: file extension for the auxiliary files (e.g. ".hdr" for .bil raster files, or ".rdc" for Idrisi .rst files)
  # ...: arguments to pass to the 'zonal' function from the 'raster' package
{
  require(raster)
  rast.files <- unzip(zip.file, list = TRUE) $ Name
  var.names <- unique(tools::file_path_sans_ext(rast.files))
  n.var <- length(var.names)
  zonal.stats <- vector("list", length(var.names))
  names(zonal.stats) <- var.names
  for (i in 1:n.var) {
    if (verbosity >= 1) message("Getting variable ", i, " of ", n.var)
    if (verbosity >= 2) message("  - unzipping file...")
    unzip(zip.file, files = paste0(var.names[i], rast.file.ext))
    if (!is.null(aux.file.ext)) unzip(zip.file, files = paste0(var.names[i], aux.file.ext))
    var.rast <- raster(paste0(var.names[i], rast.file.ext))
    if (!compareRaster(var.rast, zones.rast, stopiffalse = FALSE)) {
      if (verbosity >= 2) message("  - cropping to zones raster...")
      var.rast <- crop(var.rast, zones.rast)
    }
    if (verbosity >= 2) message("  - calculating zonal stats...")
    zonal.stats[[i]] <- raster::zonal(var.rast, zones.rast, ...)
    if (verbosity >= 2) message("  - deleting unzipped file...")
    if (delete.unzipped) {
      unlink(list.files()[grep(pattern = paste0(var.names[i], rast.file.ext), x = list.files())])
      if (!is.null(aux.file.ext)) unlink(list.files()[grep(pattern = paste0(var.names[i], aux.file.ext), x = list.files())])
    }
  }
  if (verbosity >= 1) message("Converting results to data frame...")
  zonal.stats <- as.data.frame(zonal.stats)
  zonal.stats <- subset(zonal.stats, select = c(1, seq(2, ncol(zonal.stats), 2)))
  colnames(zonal.stats)[1] <- "zone"
  colnames(zonal.stats)[-1] <- var.names
  if (verbosity >= 1) message("Finished!")
  return(zonal.stats)
}

 

Mind that you need the raster R package installed for this, and a raster map of the spatial units (zones) to which you want to extract the raster variables.Usage examples:

LGM.CCSM4.utm10 <- zonalFromZip(zip.file = "LGM/CCSM4.zip", zones.rast = utm10.rst, rast.file.ext = ".tif", aux.file.ext = NULL)
head(LGM.CCSM4.utm10)
 
WClim.utm10 <- zonalFromZip(zip.file = "bio_2-5m_bil.zip", zones.rast = utm10.rst, rast.file.ext = ".bil", aux.file.ext = ".hdr")
head(WClim.utm10)

Example for several .zip files within a folder at once (DEPRECATED – see next post’s zonalsFromZips function instead):

for (f in list.files("LGM")) {  # "LGM" is the folder containing the zip files to extract zonal stats from
  name <- paste("LGM", tools::file_path_sans_ext(f), "utm10", sep = ".")
  zonstat <- zonalFromZip(zip.file = paste("LGM", f, sep = "/"), zones.rast = u10.rst, raster.ext = ".tif", fun = "mean")
  assign(name, zonstat)
}

[presented with Pretty R]

 

Advertisements

Comment

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s