# Note: the example below requires the 'sf'-package.
if (require(sf)) {
# read vector representation of the "Mijdrecht" area
shp <- as(st_read(
dsn = system.file("maps", package = "spcosa"),
layer = "mijdrecht"), "Spatial")
# stratify into 30 strata
myStratification <- stratify(shp, nStrata = 30, nTry = 10, verbose = TRUE)
# random sampling of two sampling units per stratum
mySamplingPattern <- spsample(myStratification, n = 2)
# plot sampling pattern
plot(myStratification, mySamplingPattern)
# simulate data
# (in real world cases these data have to be obtained by field work etc.)
myData <- as(mySamplingPattern, "data.frame")
myData$observation <- rnorm(n = nrow(myData), mean = 10, sd = 1)
# design-based inference
estimate("spatial mean", myStratification, mySamplingPattern, myData["observation"])
estimate("sampling variance", myStratification, mySamplingPattern, myData["observation"])
estimate("standard error", myStratification, mySamplingPattern, myData["observation"])
estimate("spatial variance", myStratification, mySamplingPattern, myData["observation"])
estimate("scdf", myStratification, mySamplingPattern, myData["observation"])
}
Run the code above in your browser using DataLab