Apply the spatial decorrelation transformation for point-referenced data, allowing for random effects, anisotropy, partition factors, and big data methods.
decorrelate(
formula,
data,
spcov_type,
spcov_params,
xcoord,
ycoord,
algorithm = "ranger",
statistic = "RMSPE",
training,
evaluate_test,
anisotropy = FALSE,
random,
randcov_params,
partition_factor,
ordering,
local,
grid,
dense_grid,
...
)# S3 method for decorrelate_grid
tidy(x, sort_by, decreasing, ...)
A list with many elements that store information about the fitted model object:
algorithm: The machine learning algorithm used.
decorrelate_data: The output of decorrelate_data() applied to data.
fit: The fitted machine learning model object applied to the decorrelated data.
grid: If used, the grid of spatial decorrelation parameters evaluated and their corresponding
metrics when applied to the test data.
newdata: The rows of data that have NA response values and are stored as prediction data.
training: If used, the observations assigned to each training and test data set.
test: If used, a list with the lowest (absolute) mean bias (bias), mean-squared-prediction error (MSPE),
root-mean-squared-prediction error (RMSPE), and predictive R-squared (cor2).
A two-sided linear formula describing the fixed effect structure
of the model, with the response to the left of the ~ operator and
the terms on the right, separated by + operators. . on the
right-hand side represents every variable in data except the
response and the x-coordinate/y-coordinate columns (xcoord/ycoord,
or, for an sf object, the geometry column), which are never
included via . (though they may still be given explicitly).
A data frame or sf object object that contains
the variables in fixed, random, and partition_factor
as well as geographical information. If an sf object is
provided with POINT geometries, the x-coordinates and y-coordinates
are used directly. If an sf object is
provided with POLYGON geometries, the x-coordinates and y-coordinates
are taken as the centroids of each polygon.
The spatial covariance type. Available options include
"exponential", "spherical", "gaussian",
"triangular", "circular", "cubic",
"pentaspherical", "cosine", "wave",
"jbessel", "gravity", "rquad",
"magnetic", "matern", "cauchy", "pexponential",
and "none". Parameterizations of each spatial covariance type are
available in Details. Multiple spatial covariance types can be provided as
a character vector, and then decorrelate() is called iteratively for each
element and a list is returned for each model fit. The default for
spcov_type is "exponential". When spcov_type is
specified, all spatial covariance parameters are estimated.
spcov_type is ignored if spcov_params is provided.
An object from spcov_params() that contains the
spatial covariance parameters used by the spatial decorrelation transformation.
The name of the column in data representing the x-coordinate.
Can be quoted or unquoted. Not required if data is an sf object.
The name of the column in data representing the y-coordinate.
Can be quoted or unquoted. Not required if data is an sf object.
The machine learning algorithm applied. Available options
include "ranger", "randomForest", and "xgboost".
"ranger" specifies a random forest via ranger::ranger().
"randomForest" specifies a random forest via randomForest::randomForest().
"xgboost" specifies a boosted decision tree ensemble via xgboost::xgboost().
The statistic used to evaluate fit in the test data. Available options
include "bias" (mean bias), "MSPE" (mean-squared-prediction error),
"RMSPE" (root-mean-squared-prediction error)
and "cor2" (the predictive R-squared; i.e., the
squared correlation between observations and predictions).
An list controlling how the training and test data are assigned when evaluating test data performance. The following arguments detail this process:
method: The method used to evaluate test data performance.
"split" will split data up into distinct training and test sets
proportionally based on p. "cv" will split data up
via k-fold cross validation based on folds, the number of folds.
p: The proportion (a numeric vector between zero and one) of observations in data that should
be assigned to the training data. The default is 0.8, which means that
80% of the observations are assigned to the training data and 20% to the
test data. Ignored if training_index or test_index are provided.
replicate: The number of times to replicate "split" with different random training and test assignments.
folds: The number of folds to use in cross-validation. Requires method = "cv". Ignored if folds_index is specified. The default is 5, matching the 80/20 default split above.
training_index: A numeric vector that specifies which rows (i.e., indices)
of data should be assigned to the training data. If omitted, defaults
to the rows which are not already included in test_index.
test_index: A numeric vector that specifies which rows (i.e., indices)
of data should be assigned to the test data. If omitted, defaults
to the rows which are not already included in training_index.
folds_index: A numeric vector that specifies which rows
of data are associated with each cross-validation fold. Requires
method = "cv".
If omitted, training is transformed into
list(method = "split", p = 0.8, replicate = 1).
A logical indicating whether a grid should be constructed
and evaluated when spatial decorrelation parameters are known (i.e.,
spcov_params is specified, and, if random effects are included, randcov_params is specified).
If TRUE, constructs and evalutes the grid after assigning observations to
test and training data sets. If any parameters (spatial or random effects) are estimated via a grid search,
evaluate_test is set to TRUE.
A logical indicating whether (geometric) anisotropy should
be modeled. Not required if the rotate and scale parameters in spcov_params() are
0 and 1, respectively. When anisotropy is TRUE,
computational times can significantly increase. The default is FALSE.
A one-sided linear formula describing the random effect structure
of the model. Terms are specified to the right of the ~ operator.
Each term has the structure x1 + ... + xn | g1/.../gm, where x1 + ... + xn
specifies the model for the random effects and g1/.../gm is the grouping
structure. Separate terms are separated by + and must generally
be wrapped in parentheses. Random intercepts are added to each model
implicitly when at least one other variable is defined.
If a random intercept is not desired, this must be explicitly
defined (e.g., x1 + ... + xn - 1 | g1/.../gm). If only a random intercept
is desired for a grouping structure, the random intercept must be specified
as 1 | g1/.../gm. Note that g1/.../gm is shorthand for (1 | g1/.../gm).
If only random intercepts are desired and the shorthand notation is used,
parentheses can be omitted.
An object from randcov_params() that contains the
random effect variances used by the spatial decorrelation transformation.
A one-sided linear formula with a single term specifying the partition factor. The partition factor assumes observations from different levels of the partition factor are uncorrelated.
The data ordering applied. Available options
include "grts", "maxmin", "middleout",
"outsidein", "coordinate", "random", and "none".
"grts" applies ordering using a spatially balanced GRTS sample via spsurvey::grts().
"maxmin" applies maximum minimum distance ordering via GPvecchia::order_maxmin_exact().
"middleout" applies middle out ordering via GPvecchia::order_middleout().
"outsidein" applies middle out ordering via GPvecchia::order_outsidein().
"coordinate" applies middle out ordering via GPvecchia::order_coordinate(..., coordinate = c(1, 2)),
which orders from bottom-left to top-right of the spatial domain.
"random" applies a completely random ordering.
"none" applies no random ordering.
The default is "maxmin" unless there are multiple observations at a single
location, in which case the default is "grts".
A optional logical or list controlling the big data approximation.
If omitted, local is set
to TRUE or FALSE based on the sample size (the number of
non-missing observations in data) -- if the sample size exceeds 5,000,
local is set to TRUE. Otherwise it is set to FALSE.
If local is FALSE, no big data approximation
is implemented. If a list is provided, the following arguments detail the big
data approximation:
method: The big data approximation method. If method = "all",
all observations are used and size is ignored. If method = "distance",
the size data observations closest (in terms of Euclidean distance)
to the observation requiring prediction are used.
If method = "covariance", the size data observations
with the highest covariance with the observation requiring prediction are used.
If random effects and partition factors are not used in estimation and
the spatial covariance function is monotone decreasing,
"distance" and "covariance" are equivalent. The default
is "covariance".
size: The number of data observations to use when method
is "distance" or "covariance". The default is 30.
parallel: If TRUE, parallel processing via the
parallel package is automatically used. This can significantly speed
up computations even when method = "all" (i.e., no big data
approximation is used), as predictions
are spread out over multiple cores. The default is FALSE.
ncores: If parallel = TRUE, the number of cores to
parallelize over. The default is the number of available cores on your machine.
When local is a list, at least one list element must be provided to
initialize default arguments for the other list elements.
If local is TRUE, defaults for local are chosen such
that local is transformed into
list(size = 30, method = "covariance", parallel = FALSE).
An explicit grid of parameter values by which to evaluate fit. The
names of grid must contain all the names returned by decorrelate_grid(formula, data, ...).
If grid is not provided, dense_grid is a logical
which controls the density of the constructed grid to be evaluated. If
dense_grid is TRUE, a denser grid is used. If dense_grid
is FALSE, a sparser grid is used. By default, dense_grid
is FALSE when the sample size is greater than 5,000 and TRUE
otherwise.
Other arguments to the functions called by algorithm.
An object from object$grid.
Sort by a specific row in x. Fit statistics are
"bias", "MSPE", "RMSPE", and "cor2".
The default is "MSPE".
Whether sort_by should sort by decreasing order? If
sort_by = "cor2", the default is TRUE; otherwise it is FALSE.
The spatial decorrelation transformation is a preprocessing transformation that reduces the impacts of spatial dependence (i.e., covariance, correlation) on machine learning models. Predictions are made on the decorrelated scale and then recorrelated to account for spatial dependence. See Heaton et al., 2025 for details.
spcov_type Details: The correlation matrix \(R\) controls the spatial dependence structure
among observations. Parametric forms for \(R\) are given below, where \(\eta = h / range\)
for \(h\) distance between observations:
exponential: \(exp(- \eta )\)
spherical: \((1 - 1.5\eta + 0.5\eta^3) * I(h <= range)\)
gaussian: \(exp(- \eta^2 )\)
triangular: \((1 - \eta) * I(h <= range)\)
circular: \((1 - (2 / \pi) * (m * sqrt(1 - m^2) + sin^{-1}(m))) * I(h <= range), m = min(\eta, 1)\)
cubic: \((1 - 7\eta^2 + 8.75\eta^3 - 3.5\eta^5 + 0.75\eta^7) * I(h <= range)\)
pentaspherical: \((1 - 1.875\eta + 1.25\eta^3 - 0.375\eta^5) * I(h <= range)\)
cosine: \(cos(\eta)\)
wave: \(sin(\eta) / \eta * I(h > 0) + I(h = 0)\)
jbessel: \(Bj(h * range)\), Bj is Bessel-J function
gravity: \((1 + \eta^2)^{-0.5}\)
rquad: \((1 + \eta^2)^{-1}\)
magnetic: \((1 + \eta^2)^{-1.5}\)
matern: \(2^{1 - extra}/ \Gamma(extra) * \alpha^{extra} * Bk(\alpha, extra)\), \(\alpha = (2extra)^{0.5} * \eta\), Bk is Bessel-K function with order \(1/5 \le extra \le 5\)
cauchy: \((1 + \eta^2)^{-extra}\), \(extra > 0\)
pexponential: \(exp(h^{extra}/range)\), \(0 < extra \le 2\)
none: \(0\)
All spatial covariance functions are valid in one spatial dimension. All
spatial covariance functions except triangular and cosine are
valid in two dimensions. An alias for none is ie.
anisotropy Details: By default, all spatial covariance parameters except rotate
and scale as well as all random effect variance parameters
are assumed unknown, requiring estimation. If either rotate or scale
are given initial values other than 0 and 1 (respectively)
in spcov_params(), anisotropy is implicitly set to TRUE.
(Geometric) Anisotropy is modeled by transforming a covariance function that
decays differently in different directions to one that decays equally in all
directions via rotation and scaling of the original coordinates. The rotation is
controlled by the rotate parameter in \([0, \pi]\) radians. The scaling
is controlled by the scale parameter in \([0, 1]\). The anisotropy
correction involves first a rotation of the coordinates clockwise by rotate and then a
scaling of the coordinates' minor axis by the reciprocal of scale. The spatial
covariance is then computed using these transformed coordinates.
random Details: If random effects are used (the estimation method must be "reml" or
"ml"), the model
can be written as \(y = X \beta + Z1u1 + ... Zjuj + \tau + \epsilon\),
where each Z is a random effects design matrix and each u is a random effect.
partition_factor Details: The partition factor can be represented in matrix form as \(P\), where
elements of \(P\) equal one for observations in the same level of the partition
factor and zero otherwise. The covariance matrix involving only the
spatial and random effects components is then multiplied element-wise
(Hadmard product) by \(P\), yielding the final covariance matrix.
training Details: When replicate or the number of cross validation folds is at
least two, there are separate grids evaluated for each replication (or fold). Statistics in each grid are
averaged across replications (or folds) to determine a final grid ranked by statistic.
local Details: The big data approximation works by leveraging the
conditional nature of the spatial decorrelation transformation via the
Vecchia approximation. The Vecchia approximation enables efficient computation
of the conditional distribution by considering only the size most relevant
observations in the ordering (rather than using all the observations).
Observations with NA response values are removed for model
fitting, but their values can be predicted afterwards by running
predict(object).
Matthew J. Heaton, Andrew Millane, and Jake S. Rhodes. 2025. A Scalable Spatial Decorrelation Preprocessing Approach for Machine and Deep Learning. Journal of Data Science. 1-15, DOI 10.6339/25-JDS1210
decorr <- decorrelate(log_cond ~ temp, data = lake, spcov_type = "exponential")
tidy(decorr$grid)
Run the code above in your browser using DataLab