intamap-package: A package providing methods for automatic interpolation: pre-processing, parameter estimation, spatial prediction and post processing

Description

This package provides functionality for automatic interpolation of spatial data. The package was originally developed as the computational back-end of the intamap web service, but is now a stand-alone package as maintenance of the web service has ended.

Arguments

General setup

The normal work flow for working with the intamap package can best be illustrated with the following R-script. The procedure starts with reading data and meta data, then setting up an object which is used in the following functions: preprocess data, estimate parameters, compute spatial predictions, and post process them (i.e., write them out):


library(intamap)
library(sf)
# set up intamap object, either manually:
data(meuse)
coordinates(meuse) = ~x + y
data(meuse.grid)
gridded(meuse.grid) = ~x + y
obj = list(
        observations = meuse,
        predictionLocations = meuse.grid,
        targetCRS = "+init=epsg:3035",
        params = getIntamapParams()
)
class(obj) = c("idw")
# or using createIntamapObject
obj = createIntamapObject(
        observations = meuse,
        predictionLocations = meuse.grid,
        targetCRS = "+init=epsg:3035",class = c("idw")
)

# run test:
checkSetup(obj)
# do interpolation steps:
obj = preProcess(obj)
obj = estimateParameters(obj) # faster
obj = spatialPredict(obj)
obj = postProcess(obj)

Our idea is that a script following this setup will allow the full statistical analysis required for the R back-end to the automatic interpolation service, and provides the means to extend the current (over-simplistic) code with the full-grown statistical analysis routines developed by INTAMAP partners. Running the package independently under R gives the user more flexibility in the utilization than what is possible through the web-interface.

Let us look into detail what the code parts do:


library(intamap)

The command library(intamap) loads the R code of the intamap package to the current R session, along with the packages required for this (sp, gstat, akima, automap, mvtnorm, evd, MASS). After the retirement of rgdal, it is recommended to use sf-objects from the sf package. All packages need to be available to the R session, which is possible after downloading them from the Comprehensive R Network Archives (CRAN) (https://cran.r-project.org)


# set up intamap object:
obj = createIntamapObject(
        observations = meuse,
        predictionLocations = meuse.grid,
        targetCRS = "+init=epsg:3051", 
        class = "idw"
)

This code sets up a list object called obj, and assigns a class (or a group of classes) to it. This list should hold anything we need in the next steps, and the bare minimum seems to be measured point data (which will be extended to polygon data) and prediction locations, and a suggestion what to do with it. Here, the data are read from a PostGIS data base running on localhost; data base connections over a network are equally simple to set up. From the data base postgis the tables eurdep.data and inspire1km.grid are read; it is assumed that these have their SRID (spatial reference identifier) set.

The suggestion what to do with these data is put in the classes, idw. This will determine which versions of preProcess, parameterEstimate etc will be used: intamap provides methods for each of the generic functions preProcess, estimateParameters, spatialPredict, postProcess. Although it would be possible to apply two classes in this case (dataType in addition to idw), as the choice of pre- and post-processing steps tend to be data-dependent, we have tried to limit the number of classes to one for most applications.

The S3 method mechanism (used here) hence requires these versions to be called preProcess.idw, estimateParameters.idw, spatialPredict.idw, and postProcess.idw (and eventually also preProcess.eurdep and preProcess.eurdep).

To see that, we get in an interactive session:


> library(intamap)
Loading required package: sp
Loading required package: gstat
Geospatial Data Abstraction Library extensions to R successfully loaded
> methods(estimateParameters)
[1] estimateParameters.automap*         estimateParameters.copula*         
[3] estimateParameters.default*         estimateParameters.idw*            
[5] estimateParameters.linearVariogram* estimateParameters.transGaussian*  
[7] estimateParameters.yamamoto*

Now if a partner provides additional methods for BayesianKriging, one could integrate them by


class(obj) = "BayesianKriging"

and provide some or all of the functions preProcess.BayesianKriging, estimateParameters.BayesianKriging, spatialPredict.BayesianKriging, and postProcess.BayesianKriging, which would be called automatically when using their generic form (preProcess etc).

It is also possible to provide a method that calls another method. Further, for each generic there is a default method. For estimateParameter and spatialPredict these print an error message and stop, for the pre- and postprocessing the default methods may be the only thing needed for the full procedure; if no preProcess.BayesianKriging is found, preProcess.default will be used when the generic (preProcess) is called.

If a method does something, then it adds its result to the object it received, and returns this object. If it doesn't do anything, then it just passes (returns) the object it received.

To make these different methods exchangable, it is needed that they can all make the same assumptions about the contents of the object that they receive when called, and that what they return complies with what the consequent procedures expect. The details about that are given in the descriptions of the respective methods, below.

Because a specific interpolation method implemented may have its peculiar characteristics, it may have to extend these prescriptions by passing more information than described below, for example information about priors from estimateParameters to spatialPredict.

The choice between methods is usually done based on the type of problem (extreme values present, computation time available etc.). The possibility for parallel processing of the prediction step is enabled for some of the main methods. To be able to take advantage of multiple CPUs on a computer, the package doParallel must be installed, additionally the parameter nclus must be set to a value larger than 1.

Input object components

observations: object of class SpatialPointsDataFrame, containing a field value that is the target variable.
predictionLocations: object extending class Spatial, containing prediction locations.
targetCRS: character; target CRS or missing
formulaString: formula string for parameter estimation and prediction functions
params: list of parameters, to be set in getIntamapParams. These parameters include:

References

Pebesma, E., Cornford, D., Dubois, G., Heuvelink, G.B.M., Hristopulos, D., Pilz, J., Stohlker, U., Morin, G., Skoien, J.O. INTAMAP: The design and implementation f an interoperable automated interpolation Web Service. Computers and Geosciences 37 (3), 2011.