Package 'autoSTK'

Automatic Spatio-Temporal Variogram Fitting

Description

Fits a spatio-temporal variogram model to a spatio-temporal data frame using automatic parameter estimation. Supports multiple optimisers and objective functions (v2.0+).

Usage

autofitVariogramST(
  stf,
  formula,
  typestv = "sumMetric",
  candidate_model = c("Ste", "Exc", "Exp", "Wav"),
  guess_nugget = NULL,
  guess_psill = NULL,
  tlags = 0:6,
  cutoff = 20000,
  width = 500,
  aniso_method = "vgm",
  type_joint = "Exp",
  prodsum_k = NULL,
  surface = FALSE,
  measurement_error = c(0, 0, 0),
  cores = 1L,
  verbose = FALSE,
  optimizer = c("lbfgsb", "grid", "sa", "ga"),
  objective = c("WLS", "MLE"),
  n_restart = 1L,
  optimizer_control = list()
)

Arguments

stf	A spatio-temporal data frame (STFDF, STSDF, STIDF, or sftime).
formula	A formula specifying the response and explanatory variables.
typestv	Character. ST variogram model type. One of "sumMetric", "separable", "productSum", "productSumOld", "simpleSumMetric", "metric".
candidate_model	Character vector of candidate spatial/temporal variogram model families (e.g. c("Ste", "Exc", "Exp", "Wav")).
guess_nugget	Optional initial guess for the nugget. Estimated if NULL.
guess_psill	Optional initial guess for the partial sill. Estimated if NULL.
tlags	Integer vector of temporal lags.
cutoff	Numeric. Maximum spatial lag distance.
width	Numeric. Width of spatial lag bins.
aniso_method	Character. Method for estimating ST anisotropy ratio ("vgm" or "linear").
type_joint	Character. Model family for the joint variogram component.
prodsum_k	Numeric. k parameter for productSum models.
surface	Logical. If TRUE, also compute and return the variogram surface.
measurement_error	Numeric vector of length 3: spatial, temporal, joint.
cores	Integer. Number of cores for variogramST.
verbose	Logical. If TRUE, print diagnostic messages.
optimizer	Character. Optimisation strategy: "lbfgsb" (default - multi-start L-BFGS-B), "grid" (LHS grid search + L-BFGS-B refinement), "sa" (simulated annealing via optim(method = "SANN")), or "ga" (genetic algorithm via the GA package).
objective	Character. Fitting criterion: "WLS" (weighted least squares, default) or "MLE" (Gaussian log-likelihood; only for n < 500).
n_restart	Integer. Number of optimisation restarts for optimizer = "lbfgsb". Default 1 matches v1.x behaviour.
optimizer_control	Named list of extra arguments forwarded to the optimiser. Recognised keys differ by optimiser: "lbfgsb" maxit (default 2500). "grid" n_coarse (default 50), n_refine (default 30), maxit (default 2500). "sa" maxit (default 5000), temp (initial temperature, default 10), tmax (steps before cooling, default 10). "ga" popSize (default 50), maxiter (default 200), run (early-stopping generations without improvement, default 30), seed (default 42).

Value

A STVariogramFit object (list) with elements:

jointSTV

Fitted ST variogram model.

empSTV

Empirical ST variogram.

SpV

Fitted spatial marginal variogram.

TV

Fitted temporal marginal variogram.

STVsurface

(Optional) Variogram surface if surface = TRUE.

optimizer

Character. Optimiser used.

objective

Character. Objective function used.

loglik

Numeric. Log-likelihood (only when objective = "MLE").

n_obs

Integer. Number of observations (for AIC/BIC).

See Also

fitVariogramST, selectModelST, test_separability

Examples

library(spacetime)
library(sp)
library(sftime)
data(air)
stations <- st_as_sf(stations)
stations <- st_transform(stations, 'EPSG:3857')
airdf <- data.frame(PM10 = as.vector(air))
stations_full <- do.call(c, rep(stations, length(dates)))
dates_full <- rep(dates, each = nrow(stations))
rural <- cbind(airdf, time = dates_full, stations_full)
rural <- st_as_sftime(rural, sf_column_name = 'geometry')
rr <- rural[match(rural$time, dates[3001:3060], nomatch = FALSE) > 0, ]
rr <- as(as(as(rr, "STIDF"), "STFDF"), "STSDF")
rrstv <-
  autofitVariogramST(
    stf = rr,
    formula = PM10 ~ 1,
    surface = TRUE,
    cutoff = 2e6,
    width = 2e5,
    tlags = 0:6)
rrstv

---

Automatic Spatio-Temporal Kriging

Description

Convenience wrapper that fits the variogram and predicts in one call. Internally calls fitVariogramST then predictKrigeST.

Usage

autoKrigeST(
  formula,
  input_data,
  new_data,
  type_stv = "sumMetric",
  block = 0,
  model = c("Sph", "Exp", "Gau", "Ste"),
  kappa = c(0.05, seq(0.2, 2, 0.1), 5, 10),
  fix.values = c(NA, NA, NA),
  newdata_mode = "rect",
  newdata_npoints = 3000,
  GLS.model = NA,
  tlags = 0:6,
  cutoff = 20000,
  width = 500,
  forward = 6,
  predict_chunk = NULL,
  nmax = Inf,
  aniso_method = "vgm",
  type_joint = "Exp",
  prodsum_k = 0.25,
  surface = FALSE,
  start_vals = c(NA, NA, NA),
  miscFitOptions = list(),
  measurement_error = c(0, 0, 0),
  cores = 1L,
  verbose = TRUE,
  optimizer = "lbfgsb",
  objective = "WLS",
  n_restart = 1L,
  optimizer_control = list()
)

Arguments

formula	A formula (e.g. PM10 ~ 1).
input_data	Training data (STFDF, STSDF, STIDF, or sftime).
new_data	Prediction grid. If missing, generated automatically via create_new_data.ST.
type_stv	Character. ST variogram model type. Passed as typestv to autofitVariogramST.
block	Numeric. Block size for block Kriging.
model	Character vector. Candidate marginal variogram model families.
kappa	Numeric vector. Kappa values for Matern-family models.
fix.values	Numeric vector (3). Fixed initial values for nugget, range, sill.
newdata_mode	Character. "rect" or "chull".
newdata_npoints	Integer. Number of prediction grid points.
GLS.model	Variogram model for GLS sample variogram (rarely needed).
tlags	Integer vector. Temporal lags.
cutoff	Numeric. Maximum spatial lag.
width	Numeric. Spatial lag bin width.
forward	Integer. Number of future time steps to predict.
predict_chunk	Integer or NULL. Chunk size for large prediction grids.
nmax	Numeric. Maximum number of ST neighbours.
aniso_method	Character. Anisotropy estimation method.
type_joint	Character. Joint variogram component model family.
prodsum_k	Numeric. k parameter for productSum models.
surface	Logical. Return variogram surface?
start_vals	Numeric vector (3). Not used in v2.0 (kept for API compat).
miscFitOptions	List. Not used in v2.0 (kept for API compat).
measurement_error	Numeric vector (3). Measurement error components.
cores	Integer. Cores for variogramST computation.
verbose	Logical. Print progress messages.
optimizer	Character. One of "lbfgsb" (default), "grid", "sa" (simulated annealing), or "ga" (genetic algorithm).
objective	Character. "WLS" or "MLE".
n_restart	Integer. Number of optimisation restarts.
optimizer_control	List. Extra control arguments for the optimiser.

Value

An autoKrigeST object.

Examples

library(spacetime)
library(gstat)
data(air)
deair <- STFDF(stations, dates, data.frame(PM10 = as.vector(air)))
deair_rs <- deair[, 3751:3800]
## Not run:
## Not run: 
akst <- autoKrigeST(formula = PM10 ~ 1, input_data = deair_rs,
                    cutoff = 300000, width = 30000, tlags = 0:7, cores = 4)

## End(Not run)

---

Cross-validation of spatiotemporal Kriging

Description

Cross-validation of spatiotemporal Kriging

Usage

autoKrigeST.cv(
  data,
  fold_dim = c("spatial", "temporal", "random", "spacetime"),
  nfold = 10L,
  formula,
  type_stv = "sumMetric",
  block = 0,
  model = c("Sph", "Exp", "Gau", "Ste"),
  kappa = c(0.05, seq(0.2, 2, 0.1), 5, 10),
  fix.values = c(NA, NA, NA),
  tlags = 0:6,
  cutoff = 20000,
  width = 500,
  nmax = Inf,
  aniso_method = "vgm",
  type_joint = "Exp",
  prodsum_k = 0.25,
  surface = FALSE,
  start_vals = c(NA, NA, NA),
  miscFitOptions = list(),
  measurement_error = c(0, 0, 0),
  cores = 1,
  seed = 130425L,
  variogram_from_full = FALSE,
  optimizer = "lbfgsb",
  objective = "WLS",
  n_restart = 1L,
  optimizer_control = list()
)

Arguments

data	a 'STFDF'-class object
fold_dim	character. the dimension at which you want to cross-validate (spatial, temporal, and random)
nfold	integer. the number of folds. 10 as the default
formula	formula. e.g., y~1
type_stv	character. One of 'sumMetric', 'metric', 'productSum', and 'separable'
block	numeric. passed to conduct block spatiotemporal Kriging.
model	character vector. Default is c("Sph", "Exp", "Gau", "Ste"), but users can specify the list of theoretical variograms by referring gstat::vgm.
kappa	numeric vector. Kappa values tested for Matern-family variogram models.
fix.values	numeric vector. Initial values in order of nugget, range, and sill, respectively.
tlags	integer vector (increasing, preferably to be consecutive). temporal lags.
cutoff	numeric. The maximum distance at which the sample variogram will be computed.
width	numeric. The interval at which the variogram cloud will be summarized.
nmax	integer or positive infinite. The maximum number of spatiotemporal neighbors to conduct the local spatiotemporal Kriging.
aniso_method	character. One of 'vgm', 'linear', 'range', and 'metric'. Please refer to ?gstat::estiStAni.
type_joint	character. The model form of joint spatiotemporal variogram.
prodsum_k	numeric. The parameter for the case when 'productSum' is chosen for type_stv.
surface	logical. If TRUE, also return variogram surface in fitting path.
start_vals	numeric vector (3). The initial values to optimize the spatiotemporal variogram model.
miscFitOptions	list. See ?automap::autofitVariogram.
measurement_error	numeric vector of length 3: spatial, temporal, and joint measurement error terms.
cores	integer. The number of threads that will be used to compute the sample spatiotemporal variogram.
seed	integer. Random seed for fold generation.
variogram_from_full	logical. If TRUE, fit one ST variogram on full data and reuse across folds.
optimizer	character. Optimization method passed to ST variogram fitting.
objective	character. Objective for ST variogram fitting ('WLS' or 'MLE').
n_restart	integer. Number of restarts for optimizer = 'lbfgsb'.
optimizer_control	list. Extra control options forwarded to the selected optimizer.

Value

The cross-validated spatiotemporal Kriging results.

Examples

library(sp)
library(gstat)
library(spacetime)
library(stars)
library(dplyr)
data(air)
deair <- STFDF(stations, dates, data.frame(PM10 = as.vector(air)))
deair_sf <- st_as_stars(deair, crs = "+proj=longlat +ellps=sphere")
deair_sf <- st_transform(deair_sf, 3857)
deair_r <- as(deair_sf, "STFDF")
deair_r@sp@proj4string <- CRS("EPSG:3857")
deair_rs <- deair_r[, 3751:3800]
## autoKrigeST.cv test
## Not run: 
akst_cv_t <- autoKrigeST.cv(
  formula = PM10 ~ 1, data = deair_rs, nfold = 3, fold_dim = "temporal",
  cutoff = 300000, width = 30000, tlags = 0:7, cores = 8
)
akst_cv_s <- autoKrigeST.cv(
  formula = PM10 ~ 1, data = deair_rs, nfold = 3, fold_dim = "spatial",
  cutoff = 300000, width = 30000, tlags = 0:7, cores = 8
)
# akst_cv_r = autoKrigeST.cv(formula = PM10~1, data = deair_rs,  nfold = 3, fold_dim = 'random',
#                          cutoff = 300000, width = 30000, tlags = 0:7, cores = 8)
akst_cv_spt <- autoKrigeST.cv(
  formula = PM10 ~ 1, data = deair_rs, nfold = 4, fold_dim = "spacetime",
  cutoff = 300000, width = 30000, tlags = 0:7, cores = 8
)

## End(Not run)

---

Create new spatial data for interpolation

Description

Create new spatial data for interpolation

Usage

create_new_data(obj, gen_mode = "chull", npoints = 10000, return_class = "sf")

Arguments

obj	A Spatial*DataFrame.
gen_mode	character. One of 'rect' (rectangular) and 'chull' (convex hull).
npoints	integer. the number of points that will be generated
return_class	character(1). One of 'sf' and 'sp'

Value

A sf (return_class == 'sf') or SpatialPointsDataFrame (return_class == 'sp').

---

Generate a new spatiotemporal points for the spatiotemporal prediction and interpolation

Description

Generate a new spatiotemporal points for the spatiotemporal prediction and interpolation

Usage

create_new_data.ST(obj, form, gen_mode = "chull", npoints = 10000, forward = 6)

Arguments

obj	a sftime object.
form	formula.
gen_mode	character. The type of shape by the surface is generated. One of 'rect' (rectangular) and 'chull' (convex hull).
npoints	integer. the number of points that will be generated
forward	integer. the time length of the data will generate ahead of the last time point of the input data. If NULL is passed, the spatiotemporal interpolation mode in obj will be conducted.

Value

A sftime object.

---

Generate a new spatiotemporal points for the spatiotemporal prediction and interpolation (legacy sp)

Description

Generate a new spatiotemporal points for the spatiotemporal prediction and interpolation (legacy sp)

Usage

create_new_data.ST.legacy(
  obj,
  form,
  gen_mode = "chull",
  npoints = 10000,
  forward = 6
)

Arguments

obj	a ST*DF object.
form	formula.
gen_mode	character. One of 'rect' (rectangular) and 'chull' (convex hull).
npoints	integer. the number of points that will be generated
forward	integer. the time length of the data will generate ahead of the last time point of the input data. If NULL is passed, the spatiotemporal interpolation mode in obj will be conducted.

Value

A STFDF object.

---

Convert a cubble object to an sftime object

Description

Handles both nested and temporal faces of a cubble object. The resulting sftime has one row per (location × time) observation with the spatial geometry attached to every row.

Usage

cubble_to_sftime(x, time_col = NULL, key_col = NULL)

Arguments

x	A cubble_df object (from the cubble package).
time_col	Character. Name of the time column in the temporal face. If NULL (default), the function tries to detect it automatically by looking for POSIXct/Date columns.
key_col	Character. Name of the site identifier column shared by the spatial and temporal faces. Detected automatically from cubble::key_vars() when NULL.

Value

An sftime object with a geometry column (sfc_POINT).

---

Autodetect the temporal unit

Description

Autodetect the temporal unit

Usage

detect_temporal_unit(temporal)

Arguments

temporal

xts object.

Value

A character that indicates the temporal unit of the input xts object.

---

Autodetect the temporal unit in a lubridate/POSIXct object

Description

Autodetect the temporal unit in a lubridate/POSIXct object

Usage

detect_temporal_unit.generic(temporal)

Arguments

temporal

a lubridate/POSIXct object.

Value

A character that indicates the temporal unit of the input xts object.

---

Autodetect the temporal unit in a xts object

Description

Autodetect the temporal unit in a xts object

Usage

detect_temporal_unit.xts(temporal)

Arguments

temporal

a xts object.

Value

A character that indicates the temporal unit of the input xts object.

---

Estimate Initial ST Variogram Parameters from Empirical Data

Description

Derives nugget, partial sill, spatial range, and temporal range directly from the empirical spatio-temporal variogram instead of relying on fixed arithmetic heuristics. The approach is:

Usage

estimate_initial_params(
  stva_emp,
  nugget_quantile = 0.25,
  sill_quantile = 0.95,
  practical_range_target = 0.95
)

Arguments

stva_emp	An empirical StVariogram (data.frame with columns spacelag, timelag, gamma, np).
nugget_quantile	Numeric in (0, 1). Quantile of near-zero-lag gamma used as the nugget estimate. Default 0.25 (lower quartile of the first decile of spatial bins).
sill_quantile	Numeric in (0, 1). Quantile of the spatial-marginal gamma values used as the plateau estimate. Default 0.95.
practical_range_target	Numeric in (0, 1). Fraction of the sill at which the practical range is read off. Default 0.95.

Details

1. Nugget — median $\gamma$ over the first decile of spatial lags (extrapolation toward zero-lag).

2. Sill — 95th percentile of $\gamma$ values over all spatial lags at the minimum time lag (spatial marginal plateau).

3. Spatial practical range — smallest spatial lag at which the spatial marginal exceeds 95% of the estimated sill; falls back to 50% of the maximum spatial lag when the variogram never reaches the sill within the sampling window.

4. Temporal practical range — same rule applied to the temporal marginal at the minimum spatial lag.

Value

A named list:

nugget

Estimated nugget (>= 0).

psill

Estimated partial sill (sill - nugget, >= 0).

sill

Estimated total sill.

sp_range

Estimated spatial practical range.

ts_range

Estimated temporal practical range.

max_gamma

Maximum observed gamma.

max_splag

Maximum spatial lag in the empirical variogram.

max_tlag

Maximum temporal lag (numeric).

---

Fit a Spatio-Temporal Variogram

Description

A named wrapper around autofitVariogramST that makes the fit/predict separation explicit. Returns a STVariogramFit object which can be passed directly to predictKrigeST.

Usage

fitVariogramST(formula, data, ...)

Arguments

formula	A formula (e.g. PM10 ~ 1).
data	A spatio-temporal data object (STFDF, STSDF, STIDF, or sftime).
...	Arguments forwarded to autofitVariogramST.

Value

A STVariogramFit object.

See Also

predictKrigeST, autoKrigeST

---

Compute AIC and BIC for a fitted STVariogramFit (MLE-based)

Description

Compute AIC and BIC for a fitted STVariogramFit (MLE-based)

Usage

ic_stv(fit_obj, n_obs)

Arguments

fit_obj	A STVariogramFit object with an objective = "MLE" fit (i.e., fit_obj$loglik must not be NULL).
n_obs	Integer. Number of observations used for fitting.

Value

A named list: AIC, BIC, k (number of free parameters), loglik.

---

Compute the marginal spatial or temporal sample variogram

Description

Compute the marginal spatial or temporal sample variogram

Usage

marginal.variogramST(stv, bound, spatial = TRUE)

Arguments

stv	A StVariogram and data.frame object.
bound	numeric. The maximum distance that will be used to compute a spatial variogram.
spatial	boolean. if TRUE, the spatial marginal variogram will be obtained, temporal otherwise.

Value

A gstatVariogram object.

---

Sensitivity Plot of the Spatio-Temporal Anisotropy Ratio

Description

Evaluates gstat::estiStAni at a sequence of spatial intervals and plots how the estimated anisotropy ratio varies. Useful for diagnosing whether the single-interval estimate used in autofitVariogramST is stable.

Usage

plot_aniso_sensitivity(
  stva_emp,
  spatial_vgm,
  temporal_vgm,
  n_intervals = 20L,
  plot = TRUE
)

Arguments

stva_emp	Empirical ST variogram from setSTI.
spatial_vgm	Fitted spatial marginal variogram model.
temporal_vgm	Fitted temporal marginal variogram model.
n_intervals	Integer. Number of evaluation intervals.
plot	Logical. Whether to produce a base-R plot.

Value

Invisibly, a data.frame with columns interval and stAni.

---

Plot the automatically fitted variogram

Description

Plot the automatically fitted variogram

Usage

## S3 method for class 'autofitVariogram'
plot(
  x,
  plotit = TRUE,
  title = "Experimental variogram and fitted variogram model",
  ...
)

Arguments

x	A result object of autofitVariogram.
plotit	boolean. Print graph or not.
title	character. the title of the plot.
...	passed to xyplot

Value

A lattice::xyplot object.

---

Predict Using a Fitted Spatio-Temporal Kriging Model

Description

Takes a STVariogramFit object (from fitVariogramST or autofitVariogramST) and performs spatiotemporal Kriging at the locations in newdata.

Usage

predictKrigeST(
  fit,
  data,
  newdata,
  formula,
  nmax = Inf,
  predict_chunk = NULL,
  ...
)

Arguments

fit	A STVariogramFit object.
data	The original data used to fit the variogram (STFDF, STSDF, STIDF, or sftime).
newdata	Prediction grid (STFDF or sftime).
formula	Formula used during fitting.
nmax	Passed to gstat::krigeST.
predict_chunk	Optional integer; splits newdata into chunks for large grids to avoid memory issues.
...	Additional arguments forwarded to gstat::krigeST.

Value

An autoKrigeST object (list) with elements krige_output (STFDF) and var_model.

See Also

fitVariogramST, autoKrigeST

---

Automated Spatio-Temporal Covariance Model Selection

Description

Fits multiple ST covariance model types and ranks them by a user-chosen criterion.

Usage

selectModelST(
  stf,
  formula,
  candidates = c("separable", "productSum", "sumMetric", "simpleSumMetric", "metric"),
  criterion = c("WLS", "AIC", "BIC"),
  optimizer = "lbfgsb",
  cores = 1L,
  ...
)

Arguments

stf	A spatio-temporal data object (STFDF, STSDF, STIDF, or sftime).
formula	A formula (e.g., PM10 ~ 1).
candidates	Character vector of ST model types to try. Defaults to all six supported types.
criterion	One of "AIC", "BIC" (require objective = "MLE"), "CV_RMSE" (not yet implemented in v2.0; use "WLS" as a fast proxy), or "WLS".
optimizer	Passed to autofitVariogramST.
cores	Passed to autofitVariogramST.
...	Additional arguments forwarded to autofitVariogramST.

Value

An STModelSelection object (list) with elements best_model, comparison (data.frame), and all_fits.

---

A convenience function for the sample spatiotemporal variogram

Description

A convenience function for the sample spatiotemporal variogram

Usage

setSTI(
  stf,
  formula,
  tlags = 0:6,
  cutoff = 30000,
  width = 1000,
  assumeRegular = FALSE,
  pseudo = TRUE,
  na.omit = TRUE,
  wireframe = FALSE,
  plot3d = FALSE,
  cores = 1
)

Arguments

stf	A ST*DF object.
formula	formula (inherits the same parameter in variogramST)
tlags	temporal lags to compute semivariance (inherits the same parameter in variogramST)
cutoff	the maximum bound of the set of spatial lags (inherits the same parameter in variogramST)
width	integer (1). spatial lag (inherits the same parameter in variogramST)
assumeRegular	Boolean. Assuming regular grid?
pseudo	Boolean. See ?gstat::variogramST
na.omit	Boolean. Omit NA values.
wireframe	Boolean. Whether you plot a StVariogram in wireframe or not. If not, the return will be in class of data.frame, not a list
plot3d	Boolean. Wheter you make a three-dimensional graph with rgl package
cores	Integer. Number of threads passed to gstat::variogramST.

Value

Depends on the arguments wireframe (if TRUE, list of length 2) and plot3d (if TRUE, list of length 3), a StVariogram object otherwise.

---

Compute parameter bounds for vgmST model optimisation

Description

Compute parameter bounds for vgmST model optimisation

Usage

st_param_bounds(
  model_template,
  stva_emp,
  sill_scale = 2,
  range_scale = 3,
  ani_scale = 20
)

Arguments

model_template	A vgmST object (initial guess).
stva_emp	Empirical ST variogram (StVariogram/data.frame) from setSTI or gstat::variogramST.
sill_scale	Numeric (2). Upper bound multiplier applied to max(stva_emp$gamma) for sill/psill/nugget parameters.
range_scale	Numeric (3). Upper bound multiplier applied to the maximum spatial lag for range parameters.
ani_scale	Numeric (20). Controls the breadth of the stAni search interval relative to the data extent.

Value

A named list with elements lower and upper (both named numeric vectors matching gstat::extractPar(model_template)).

---

Likelihood Ratio Test for Spatiotemporal Separability

Description

Tests whether a separable covariance structure is adequate versus the more general sumMetric model using a likelihood ratio test (LRT). Both models are fitted with objective = "MLE".

Usage

test_separability(stf, formula, ...)

Arguments

stf	Spatio-temporal data (STFDF, STSDF, STIDF, or sftime).
formula	Formula, e.g. PM10 ~ 1.
...	Additional arguments forwarded to autofitVariogramST (e.g. tlags, cutoff, width, cores).

Value

A list with elements:

statistic

LRT chi-squared statistic.

df

Degrees of freedom (difference in number of parameters).

p_value

p-value from chi-squared distribution.

separable_fit

Fitted STVariogramFit for the separable model.

summetric_fit

Fitted STVariogramFit for the sumMetric model.

---