RSAGA is a package providing access to geoprocessing and terrain analysis functions of SAGA-GIS (System for Automated Geospatial Analysis; Conrad et al., 2015)1, and R functions for spatial data manipulation and applying functions to stacks of grids. This coupling supports statistical geocomputing by which terrain analysis may be executed in the familiar R environment. This tutorial provides an introduction to RSAGA, with examples to guide initializing the geoprocessing environment, accessing SAGA and ASCII spatial data, function calls between RSAGA and SAGA, and RSAGA tools to apply model predictions to grids.
This vignette illustrates the use of RSAGA for terrain analysis following an application modeling landslide susceptibility in the tropical Andes (Muenchow et al., 2012). The presence or absence of landslide initiation at 1535 sampling points in the natural parts of the perhumid Reserva Biologica San Francisco (RBSF) area are entered as the response to terrain variables in a generalized additive model (GAM). From a provided digital elevation model (DEM; figure 1) (Jordan et al., 2005; Ungerechts, L. 2010)2 , RSAGA is used to calculate the relevant terrain indices, transform grid data, and extract values to a data frame for statistical modeling. The GAM is demonstrated using package gam (Hastie, 2015), and model predictions applied to the collection of terrain grids using RSAGA grid tools.
RSAGA can be considered a collection of geocomputing functions in three general components:
RSAGA core provides access to the SAGA
command-line interpreter, converting R commands to system calls to use
SAGA. The rsaga.geoprocessor
function allows the user to
specify a list of arguments to be provided to any command-line SAGA
module. Navigating SAGA help and other module information can be
accomplished using rsaga.get.modules
and
rsaga.get.usage
functions, which print the console output
from various SAGA command-line calls.
RSAGA modules are a more R-convenient interface of commonly-used SAGA modules, written as R functions. These use the geoprocessor to create SAGA command-line calls, with an R interface including help and usage.
RSAGA grid tools provides additional access to geocomputing and data manipulation. These include moving window and local grid functions, data import and export, and “picking” functions to extract grid and shapefile data to an R data frame.
Running the geoprocessing functions requires that SAGA is installed
and can be found by RSAGA. The rsaga.env()
function creates
a list with system-dependent information on SAGA path, module path and
data directory. The function should be called with the path to the SAGA
executable:
<- rsaga.env() env
The function tries to find SAGA in several default paths if no path parameter is provided.
The object env
may now be used to direct RSAGA functions
to use this specific version of SAGA. Version 8.4.1
will be
used for this tutorial by specifying env = env
in functions
that use this SAGA environment.
RSAGA interfaces with the command line interpreter of SAGA GIS. In
most cases, rsaga.*
functions take specified arguments and
translate these to a system call that uses SAGA directly. For instance,
slope may be calculated from a DEM by the SAGA command line call:
saga_cmd ta_morphometry "Slope, Aspect, Curvature" -ELEVATION "C:\WORKSPACE\dem.sgrd"
-SLOPE "C:\WORKSPACE\slope.sgrd"
This is accomplished in RSAGA by running the geoprocessor function:
rsaga.geoprocessor(lib = "ta_morphometry", module = "Slope, Aspect, Curvature",
param = list(ELEVATION = paste(getwd(),"/dem.sgrd", sep = ""),
SLOPE = paste(getwd(),"/slope.sgrd", sep = "")),
env = env)
where
lib
is the library of SAGA modules,module
is the function to be used,param
is the list of parameters required by the command
line interpreter,env
is a SAGA geoprocessing environment containing
information on SAGA and the path to SAGA modulesTo simplify the call, geoprocessing can also be performed through
rsaga.*
functions:
rsaga.slope(in.dem = "dem", out.slope = "slope", env = env)
These functions translate arguments to the geoprocessor function,
which in turn creates a SAGA command-line call to execute the module.
Options for the geoprocessor such as display.command
and
show.output.on.console
may be specified in these
functions.
We hope you find that SAGA modules with corresponding
rsaga.*
functions have detailed help files that can be
called, e.g.:
help(rsaga.slope)
However, you may wish to access a SAGA module that does not have an RSAGA function. You may also wish to explore changes in module parameters and deprecated modules between versions of SAGA. Available libraries and modules can be found:
# A character string of available libraries:
rsaga.get.libraries(path = env$modules)
# A list of modules in a library:
rsaga.get.modules(libs = "ta_morphometry", env = env)
The usage, or command line arguments, of a particular module can be found, printing a character string to the console that lists arguments and choices for the selected module:
rsaga.get.usage(lib = "ta_morphometry", module = "Slope, Aspect, Curvature", env = env)
# Compare module parameters between versions:
rsaga.get.usage(lib = "ta_morphometry", module = "Slope, Aspect, Curvature",
env = rsaga.env(path = "C:/SAGA-GIS/saga_2.1.0_x64"))
rsaga.get.usage(lib = "ta_morphometry", module = "Slope, Aspect, Curvature",
env = rsaga.env(path = "C:/SAGA-GIS/saga_2.2.0_x64"))
In most cases RSAGA will access SAGA .sgrd
format grids
and ASCII grids. For this example, .rda
-compatible version
of an ASCII DEM as a large list is provided in the
landslides
dataset. This object is what would be created by
the RSAGA function read.ascii.grid()
. Spatial sampling
points are included as a data frame with x
and
y
coordinates. These are loaded:
data(landslides)
To use the provided DEM with SAGA through RSAGA, it can be converted
to the .sgrd
format:
write.sgrd(data = dem, file = "dem", header = dem$header,
env = env) # write.sgrd and read.sgrd use SAGA, and should specify 'env'
Or to an ASCII .asc
grid:
write.ascii.grid(data = dem, file = "dem", header = dem$header)
It is important to note that these functions, as with most RSAGA
functions, do not store grids in the R session; rather, the processing
takes place within SAGA. Grids in various workspaces can be accessed or
written by specifying the full path in the argument of that grid. For
example, in the rsaga.slope
call above, reading the DEM
required for the SAGA calculation and writing the output slope model to
a different directory can be specified, e.g.:
rsaga.slope(in.dem = "C:/InData/dem", out.slope = "C:/OutData/slope", env = env)
If no path is specified, RSAGA will read and write grids to the current R working directory.
With the DEM now as a SAGA grid, RSAGA may then be used to process terrain and morphometric variables to be used, for example, in this landslide analysis. Slopes and curvatures are calculated by the Zevenbergen & Thorne (1987) method (figure 2):
# By individual function calls:
rsaga.slope("dem", "slope", method = "poly2zevenbergen", env = env)
rsaga.plan.curvature("dem", "cplan", method = "poly2zevenbergen", env = env)
rsaga.profile.curvature("dem", "cprof", method = "poly2zevenbergen", env = env)
# By one function that calculates each of the terrain parameters:
rsaga.slope.asp.curv("dem", out.slope = "slope",
out.cprof = "cprof", out.cplan = "cplan",
method = "poly2zevenbergen",
env = env)
Contributing area (carea
) is calculated through
rsaga.topdown.processing
using the multiple flow direction
method (Quinn et al., 1991; figure 2):
rsaga.topdown.processing("dem", out.carea = "carea", method = "mfd", env = env)
In some cases, it may be desirable to use a transformation of a
terrain variable for modeling. For example, the catchment area may be
better expressed once log-transformed. Such transformations can be
accomplished in RSAGA with rsaga.grid.calculus
. The
formula
in grid calculus is very flexible. Input grids are
represented in the formula by letters a
, b
,
etc. Any of the arithmetic operators and functions listed in
the help file may be used to calculate the out.grid
. Here,
the formula ~ log(a)
will calculate the logarithm of the
contributing area:
rsaga.grid.calculus(in.grids = "carea", out.grid = "log10_carea",
formula = ~ log(a), env = env)
The grid data may be extracted spatially to a data frame coincident
with sampling points. The raster values may be “picked” by
pick.from.*
functions (technically interpolated by the
nearest neighbour method) and added to the landslides
data
frame:
# Pick grid values; add to landslides data.frame:
<- pick.from.saga.grid(landslides, "slope", varname = "slope", env = env,
landslides X.name = "x", Y.name = "y")
<- pick.from.saga.grid(landslides, "cprof", varname = "cprof", env = env,
landslides X.name = "x", Y.name = "y")
<- pick.from.saga.grid(landslides, "cplan", varname = "cplan", env = env,
landslides X.name = "x", Y.name = "y")
<- pick.from.saga.grid(landslides, "carea", varname = "carea", env = env,
landslides X.name = "x", Y.name = "y")