RStoolbox
is an R package providing a wide range of
tools for your every-day remote sensing processing needs. The available
tool-set covers many aspects from data import, pre-processing, data
analysis, image classification and graphical display.
RStoolbox
builds upon the terra
package, which
makes it suitable for processing large data-sets even on smaller
workstations.
Find out more on the RStoolbox
webpage.
The package is available on CRAN and can be installed as usual via
install.packages("RStoolbox")
To install the latest version from GitHub you need to have r-base-dev (Linux) or Rtools (Windows) installed. Then run the following lines:
library(devtools)
install_github("bleutner/RStoolbox")
RStoolbox
implements a variety of remote
sensing methods and workflows. Below are a few examples to get
started. Further examples can be found in the documentation
of the respective functions.
The example below shows an unsupervised classification workflow based on kmeans clustering:
library(RStoolbox)
# unsupervised classification with 3 classes
<- unsuperClass(lsat, nClasses = 3)
uc
# plot result using ggplot
ggR(uc$map, geom_raster = T, forceCat = T) +
scale_fill_manual(values = c("darkgreen", "blue", "sandybrown"))
If training data are available, e.g. labeled polygons,
RStoolbox
can be used to conduct a supervised
classification. The workflow below employs randomForest to
train a classification model:
library(RStoolbox)
library(caret)
library(randomForest)
library(ggplot2)
library(terra)
# example: training data from digitized polygons
<- readRDS(system.file("external/trainingPolygons_lsat.rds", package="RStoolbox"))
train
# plot input data
ggRGB(lsat, r = 3, g = 2, b=1, stretch = "lin") +
geom_sf(data = train, aes(fill = class)) +
scale_fill_manual(values = c("yellow", "sandybrown", "darkgreen", "blue"))
#> Coordinate system already present. Adding new coordinate system, which will
#> replace the existing one.
# fit random forest (splitting training into 70\% training data, 30\% validation data)
<- superClass(lsat, trainData = train, responseCol = "class",
sc model = "rf", tuneLength = 1, trainPartition = 0.7)
# print model performance and confusion matrix
$modelFit
sc#> [[1]]
#> TrainAccuracy TrainKappa method
#> 1 0.9992293 0.9988032 rf
#>
#> [[2]]
#> Cross-Validated (5 fold) Confusion Matrix
#>
#> (entries are average cell counts across resamples)
#>
#> Reference
#> Prediction cleared fallen_dry forest water
#> cleared 141.6 0.0 0.0 0.0
#> fallen_dry 0.0 22.0 0.0 0.0
#> forest 0.4 0.0 255.0 0.0
#> water 0.0 0.0 0.0 99.4
#>
#> Accuracy (average) : 0.9992
# plotting: convert class IDs to class labels (factorize) and plot
<- as.factor(sc$map)
r levels(r) <- data.frame(ID = 1:4, class_supervised = levels(train$class))
ggR(r, geom_raster = T, forceCat = T) + scale_fill_manual(values = c("yellow", "sandybrown", "darkgreen", "blue"))
Created on 2024-04-19 with reprex v2.1.0
RStoolbox
offers spectral unmixing by implementing the
Multiple Endmember Spectral Mixture Analysis (MESMA) approach for
estimating fractions of spectral classes, such as spectra of surfaces or
materials, on a sub-pixel scale. The following workflow shows a simple
Spectral Mixture Analysis (SMA) with single endmembers per class,
extracted from the lsat
example image by cell id:
library(RStoolbox)
library(terra)
# to perform a SMA, use a single endmember per class, row by row:
<- data.frame(lsat[c(5294, 47916)])
em rownames(em) <- c("forest", "water")
# umix the lsat image
<- mesma(img = lsat, em = em)
probs plot(probs)
Instead, one can define multiple endmembers per class to conduct a Multiple Endmember Spectral Mixture Analysis (MESMA):
library(RStoolbox)
library(terra)
# to perform a MESMA, use multiple endmembers per class, differntiating them
# by a column named 'class':
<- rbind(
em data.frame(lsat[c(4155, 17018, 53134, 69487, 83704)], class = "forest"),
data.frame(lsat[c(22742, 25946, 38617, 59632, 67313)], class = "water")
)
# unmix the lsat image
<- mesma(img = lsat, em = em)
probs plot(probs)
# MESMA can also be performed on more than two endmember classes:
<- rbind(
em data.frame(lsat[c(4155, 17018, 53134, 69487, 83704)], class = "forest"),
data.frame(lsat[c(22742, 25946, 38617, 59632, 67313)], class = "water"),
data.frame(lsat[c(4330, 1762, 1278, 1357, 17414)], class = "shortgrown")
)
# unmix the lsat image
<- mesma(img = lsat, em = em)
probs plot(probs)
RStoolbox
comes with a suite of pre-processing
functions, including cloudMask
to identify clouds in
optical satellite imagery:
library(ggplot2)
# lsat example scene, with two tiny clouds in the east
ggRGB(lsat, stretch = "lin")
# calculate cloud index
<- cloudMask(lsat, blue = 1, tir = 6)
cldmsk ggR(cldmsk, 2, geom_raster = TRUE)
# mask by threshold, region-growing around the core cloud pixels
<- cloudMask(cldmsk, threshold = 0.1, buffer = 5)
cldmsk_final
## plot cloudmask
ggRGB(lsat, stretch = "lin") +
ggR(cldmsk_final[[1]], ggLayer = TRUE, forceCat = TRUE, geom_raster = TRUE) +
scale_fill_manual(values = c("red"), na.value = NA)
#> Warning: Removed 88752 rows containing missing values or values outside the scale range
#> (`geom_raster()`).
With radCor
, users can compute radiometric and simple
atmospheric corrections (based on dark object substraction):
library(terra)
# import Landsat meta data
<- system.file("external/landsat/LT52240631988227CUB02_MTL.txt", package="RStoolbox")
mtlFile <- readMeta(mtlFile)
metaData <- stackMeta(mtlFile)
lsat_t
# convert DN to top of atmosphere reflectance and brightness temperature
<- radCor(lsat_t, metaData = metaData, method = "apref")
lsat_ref
# correct DN to at-surface-reflecatance with DOS (Chavez decay model)
<- radCor(lsat_t, metaData = metaData)
lsat_sref
# correct DN to at-surface-reflecatance with simple DOS and automatic haze estimation
<- estimateHaze(lsat_t, hazeBands = 1:4, darkProp = 0.01, plot = FALSE)
hazeDN <- radCor(lsat_t, metaData = metaData, method = "sdos",
lsat_sref hazeValues = hazeDN, hazeBands = 1:4)
# plot result
ggRGB(lsat_sref, r = 3, g = 2, b = 1, stretch = "lin")
Created on 2024-04-19 with reprex v2.1.0