Selection of representative observations
The selection of observations from an existing sample with rassta results in one observation per stratification unit in a set. This observation is termed the representative response observation. The representative response observation is assumed to best reflect the control that the landscape configuration represented by a stratification unit exerts on response phenomena. Currently, rassta allows to select a representative response observation for each stratification unit in a set based on the following methods:
- mls: selects the observation at the location with the maximum landscape similarity value.
- mrv: selects the observation whose response value is the median of all the values.
- random: selects an observation at random.
The first two methods represent a stratified non-probability sampling, while the third method represents a stratified random sampling. Moreover, the first method is biased towards maximizing the representativeness of landscape configurations because it uses landscape similarity layers. The codes below demonstrate the selection of representative response observations with observation(). The observations will be selected from a sample of n = 195 based on the mls method.
- Creation of multi-layer SpatRaster with landscape similarity to stratification units.
# Load rassta and terra packages
library(rassta)
library(terra)
# Single-layer SpatRaster of stratification units
su <- rast(paste(d, "/su.tif", sep = ""))
# Multi-layer SpatRaster with spatial signatures of classification units
clim.sig <- rast(list.files(d, pattern = "climate_", full.names = TRUE)) # For climatic units
mat.sig <- rast(list.files(d, pattern = "material_", full.names = TRUE)) # For material units
terr.sig <- rast(list.files(d, pattern = "terrain_", full.names = TRUE)) # For terrain units
# Landscape similarity to stratification units
su.ls <- similarity(su.rast = su, sig.rast = c(clim.sig, mat.sig, terr.sig),
su.code = list(climate = c(1, 1),
material = c(2, 2),
terrain = c(3, 3)
)
)- Selection of the representative response observation for each stratification unit. The response variable in this example is soil organic carbon (SOC).
# SpatVector with SOC observations for stratification units
socobs <- vect(paste(d, "/soc.shp", sep = ""))
# Representative SOC observation for each stratification unit
su.obs <- observation(su.rast = su, obs = socobs, col.id = 1, col.resp = 2,
method = "mls", ls.rast = su.ls$landsim
)
# Information about sample and representative SOC observations
socobs
#> class : SpatVector
#> geometry : points
#> dimensions : 195, 2 (geometries, attributes)
#> extent : 284252, 292532, 153471, 161811 (xmin, xmax, ymin, ymax)
#> source : soc.shp
#> coord. ref. : NAD83 / Washington South (EPSG:32149)
#> names : ID soc
#> type : <int> <num>
#> values : 1 6.603
#> 2 5.701
#> 3 5.594
su.obs$su_repobs.sp
#> class : SpatVector
#> geometry : points
#> dimensions : 56, 4 (geometries, attributes)
#> extent : 284252, 292532, 153651, 161811 (xmin, xmax, ymin, ymax)
#> coord. ref. : NAD83 / Washington South (EPSG:32149)
#> names : ID soc SU land_sim
#> type : <int> <num> <int> <num>
#> values : 1 6.603 111 80.67
#> 6 5.443 112 85
#> 10 4.609 114 77.33- The spatial and statistical distribution of the sample and representative response observations can be plotted as demonstrated in the code below.
# Set graphics arrangement
par(mfrow = c(1,2))
# Plot stratification units and response observations
plot(su, type = "classes", col = hcl.colors(56, "spectral"), legend = FALSE,
mar = c(3, 2, 3, 1.5), main = paste("Soil Organic Carbon Observations"),
fun = function() c(points(socobs, pch = 21, bg = rgb(0,1,0,1)),
points(su.obs$su_repobs.sp, pch = 21, bg = rgb(0,0,1,1))
)
)
# Set new graphics arrangement
par(mar = c(2, 1.5, 1.5, 1.5))
# Plot histogram of soil organic carbon values from all observations
hist(socobs$soc, 4, col = rgb(0,1,0,0.8), main = "Soil Organic Carbon (%)", xlab = "")
# Plot histogram of soil organic carbon values from representative observations
hist(su.obs$su_repobs.sp$soc, 4, add = T, col = rgb(0,0,1,0.9))
# Add legend
legend("topright", legend = c("Sample", "Representative Observations"),
col = c(rgb(0,1,0,0.8), rgb(0,0,1,0.9)), pch = 20, bty = "n", pt.cex = 2,
cex = 0.6, text.col = "black", horiz = F, inset = c(0, 0.05)
)
