intamapInteractive (version 1.2-6)

ssaOptim: Spatial simulated annealing (SSA) for optimization of sampling designs using a geostatistical measure of spatial prediction error

Description

Spatial simulated annealing uses slight perturbations of previous sampling designs and a random search technique to solve spatial optimization problems. Candidate measurement locations are iteratively moved around and optimized by minimizing the mean universal kriging variance (calculateMukv). The approach relies on a known, pre-specified model for underlying spatial variation (variogramModel).

Usage

ssaOptim(observations, predGrid, candidates, action, nDiff, model,
         nr_iterations, plotOptim = TRUE, formulaString = NULL, 
         coolingFactor = nr_iterations/10, covariates = "over", fun, ...)

Value

SpatialPointsDataFrame with optimized locations

Arguments

observations

object of class Spatial with coordinates and possible covariates

predGrid

object of class Spatial* used when method = "ssa". predGrid should contain the coordinates of prediction locations for optimization. Usually predGrid is a SpatialGrid / SpatialPixels or a SpatialGridDataFrame / SpatialPixelsDataFrame when independent covariate predictor variables are used

candidates

candidates is the study area of class SpatialPolygonsDataFrame

action

character string indicating which type of action to perform: "add" to add new measurement stations to the existing network or "del" to delete existing stations

nDiff

number of stations to add or delete

model

variogram model:object of class variogramModel, as generated by vgm

nr_iterations

number of iterations to process before stopping. The default coolingFactor is also a function of number of iterations.

plotOptim

logical; if TRUE, creates a plot of the result as optimization progresses; TRUE by default

formulaString

formula that defines the dependent variable as a linear model of independent variables; suppose the dependent variable has name z, for ordinary and simple kriging use the formula z~1; for universal kriging, suppose z is linearly dependent on x and y, use the formula z~x+y. The formulaString defaults to "value~1" if value is a part of the data set. If not, the first column of the data set is used.

coolingFactor

variable defining how fast the algorithm will cool down. With SSA, worsening designs are accepted with a decreasing probability (generally set to p $<$ 0.2 to avoid selection of local minima). The coolingFactor dictates the rate at which p decreases to zero. Commonly p is set to exponentially decrease or cool as a function of number of iterations to ensure convergence (Brus & Heuvelink 2007, Melles et al. 2011). Smaller numbers give quicker cooling; higher numbers give slower cooling.

covariates

character string defining whether possible covariates should be found by "over" or "krige", see also details below

fun

Alternative objective function for optimization, the input and output should match the ones of (calculateMukv (except for fun)

...

other arguments to be passed on to lower level functions

Author

O. Baume, S.J. Melles, J. Skoien

Details

The default version of this function applies spatial simulated annealing for optimization with the MUKV criterion (calculateMukv). With covariates, the function takes a universal kriging model into account. When optimizing a sampling design using SSA and criterion = "mukv", measurement values at new sampling locations are not required in order to calculate prediction error variance (criterion = "mukv"). This attractive property allows one to estimate the kriging prediction error variance prior to collecting the data (i.e., the dependent variable is unknown at new candidate locations), and it is this property that is used in the SSA optimization procedure after (Brus & Heuvelink 2007, Melles et al. 2011).

A stopping criterion countMax is implemented in lower level functions to end the optimization procedure after 200 search steps have occurred without an improvement in the design. If this stopping criterion is reached before nr_iterations, SSA will terminate.

method = "ssa" with criterion = "mukv" makes it possible to assume a linear relationship between independent variables in predGrid and dependent variables at observation locations using universal kriging (krige). However, a correct estimate of mean universal kriging variance requires that the independent covariate variables be known at candidate locations. Thus it is necessary to have complete spatial coverage for all covariate predictors in the linear model. Covariate information must be available at both new candidate measurement locations and prediction locations. This is not possible (except for the measurement locations themselves). Instead, these are estimated from the prediction locations.

There are two possible methods to attain information on covariates at the candidate locations, and the method can be set using the argument covariates: over and krige. over finds the value of covariates at new locations by overlaying candidate locations on the prediction grid and taking the value of the nearest neighbour. The second method uses kriging to estimate covariate values at new locations from predGrid. The first method is generally faster, the second method is most likely more exact, particularly if the resolution of predGrid is low in relation to the spatial correlation lengths of the covariates. Both methods are approximations that may influence the criterion used for optimization with increasing numbers of points added.

It is possible to submit an alternative function fun as objective function. This function should take at least the observation locations and the predGrid as input, and return a value which should be minimized. See also calculateMukv for more information about arguments to this function.

References

D. J. Brus, G. B. M. Heuvelink (2007). Optimization of sample patterns for universal kriging of environmental variables, Geoderma, 138: 86-95 (2007).

S. J. Melles, G. B. M. Heuvelink, C. J. W. Twenhofel, U. Stohlker (2011). Optimizing the spatial pattern of networks for monitoring radioactive releases, Computers and Geosciences, 37: 280-288 (2011).

Examples

# \donttest{
# load data:
library(gstat)
data(meuse)
coordinates(meuse) = ~x+y
data(meuse.grid)
coordinates(meuse.grid) = ~x+y
gridded(meuse.grid) = TRUE
predGrid = meuse.grid

# estimate variograms (OK/UK):
vfitOK = fit.variogram(variogram(zinc~1, meuse), vgm(1, "Exp", 300, 1))
vfitUK = fit.variogram(variogram(zinc~x+y, meuse), vgm(1, "Exp", 300, 1))
vfitRK = fit.variogram(variogram(zinc~dist+ffreq+soil, meuse), vgm(1, "Exp", 300, 1))

# study area of interest:
bb = bbox(predGrid)
boun = SpatialPoints(data.frame(x=c(bb[1,1],bb[1,2],bb[1,2],bb[1,1],bb[1,1]),
                                y=c(bb[2,1],bb[2,1],bb[2,2],bb[2,2],bb[2,1])))
Srl = Polygons(list(Polygon(boun)),ID = as.character(1))
candidates = SpatialPolygonsDataFrame(SpatialPolygons(list(Srl)),
                                      data = data.frame(ID=1))

# add 20 more points assuming OK model (SSA method):
optimOK <- ssaOptim(meuse, meuse.grid, candidates = candidates, covariates = "over",
            nDiff = 20, action = "add", model = vfitOK, nr_iterations = 10000, 
            formulaString = zinc~1, nmax = 40, countMax = 200)

# add 20 more points assuming UK model (SSA method):
optimUK <- ssaOptim(meuse, meuse.grid, candidates = candidates, covariates = "over",
            nDiff = 20, action = "add", model = vfitUK, nr_iterations = 10000, 
            formulaString = zinc~x+y, nmax = 40, countMax = 200)

# add 20 more points with auxiliary variables (SSA method):
optimRK <- ssaOptim(meuse, meuse.grid, candidates = candidates, covariates = "over",
            nDiff = 20, action = "add", model = vfitRK, nr_iterations = 10000, 
            formulaString = zinc~dist+ffreq+soil, nmax = 40, countMax = 200)
   # }