Simulate a spatio-temporal Gaussian random field on a regular grid

Description

Simulates a space-time Gaussian random field on a regular (x,y,t) grid. The field is returned as a 3D array and can be used as a latent field for log-Gaussian Cox process (LGCP) simulation.

Usage

Gauss.st.F(
  xlim = c(0, 1),
  ylim = c(0, 1),
  tlim = c(0, 1),
  par1 = c(1, 0.05),
  par2 = c(1, 0.06),
  sigmas = c(0.5, 0.5, 1),
  grid = c(15L, 15L, 10L)
)

Arguments

Details

The simulated field is a weighted sum of three independent Gaussian components:

Z(x,y,t) = \sigma_1 Z_s(x,y) + \sigma_2 Z_t(t) + \sigma_3 Z_{st}(x,y,t),

where Z_s is a purely spatial field, Z_t is a purely temporal field, and Z_{st} is a spatio-temporal field with separable exponential covariance in space and time.

The function uses mvrnorm for multivariate normal simulation and rdist to compute pairwise distances for covariance matrix construction.

Spatial and temporal covariances are exponential. The spatio-temporal component uses a separable covariance C_{st}((u,t),(u',t')) = C_s(u,u') C_t(t,t'). Simulation is performed via Cholesky factors without constructing the full (nx*ny*nt) \times (nx*ny*nt) covariance matrix.

Value

A list with components:

Z

Numeric array of dimension c(nx, ny, nt) containing simulated field values.

xcoord

Numeric vector of length nx with x-grid coordinates.

ycoord

Numeric vector of length ny with y-grid coordinates.

tcoord

Numeric vector of length nt with time-grid coordinates.

Author(s)

Mohammad Ghorbani mohammad.ghorbani@slu.se
Nafiseh Vafaei nafiseh.vafaei@slu.se

References

Ghorbani M., Vafaei N., Dvořák J., Myllymäki M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics and Data Analysis, 161, 107245.

See Also

mvrnorm, rdist

Examples

 if (requireNamespace("MASS", quietly = TRUE) && requireNamespace("fields", quietly = TRUE)) {
  set.seed(1)
  field <- Gauss.st.F(
    xlim = c(0, 1), ylim = c(0, 1), tlim = c(0, 1),
    par1 = c(1, 0.05), par2 = c(1, 0.06),
    sigmas = c(0.5, 0.5, 1),
    grid = c(15, 15, 10)
  )
# Inspect dimensions and visualize one time slice
dim(field$Z)
image(field$xcoord, field$ycoord, field$Z[, , 1],
      main = "Gaussian Random Field (t = 1)",
      col = RColorBrewer::brewer.pal(11, "Spectral"))
  }

Compute S-based test function for testing the null hypothesis of first-order separability

Description

Computes kernel-based estimates of the spatio-temporal intensity and related separability diagnostics, either on a regular spatio-temporal grid ("pixels") or at the observed event locations ("points").

Usage

S.based.functions(
  X,
  s.region,
  t.region,
  owin = NULL,
  at = c("pixels", "points"),
  n.grid = c(25L, 25L, 20L),
  epsilon = NULL,
  delta = NULL,
  output = "all"
)

Arguments

Details

The function is a wrapper that (i) validates inputs, (ii) computes bandwidths and Gaussian edge-correction masses via calc.bandwidths.and.edgecorr, and (iii) delegates the actual estimation to estimate.intensity.pixel or estimate.intensity.point.

When at="points", spatial/temporal profiles such as S.space and S.time are typically not defined and are returned as NULL by the pointwise routine.

Value

If output="all", returns the full list produced by the chosen computation routine, augmented with epsilon and delta.

If output is a single component name, returns a list with:

S

The requested component.

epsilon

Spatial bandwidth used.

delta

Temporal bandwidth used.

x,y,t

Coordinates returned by the underlying routine (may be NULL).

Examples

X <- cbind(stats::runif(50), stats::runif(50), stats::runif(50))
s.region <- matrix(c(0,0, 1,0, 1,1, 0,1), ncol = 2, byrow = TRUE)
t.region <- c(0, 1)
res_all <- S.based.functions(X, s.region, t.region, at = "points")
res_Sfun <- S.based.functions(X, s.region, t.region, at = "points", output = "S.fun")

Block permutation of the temporal component in a spatio-temporal point pattern

Description

Permutes the temporal component of a spatio-temporal dataset in a block-wise manner while keeping the spatial coordinates fixed. This is used to generate permuted replicates under the null model of first-order separability.

Usage

block.permut(nblocks, X, nperm = 1999)

Arguments

Details

The function first orders the events by time and partitions the ordered sequence into nblocks consecutive blocks of equal size. The block labels are permuted (excluding the identity permutation), and the time values are reassigned according to the permuted block order.

If nrow(X) is not divisible by nblocks, the last nrow(X) %% nblocks events are not included in the block permutation and are appended unchanged to each permuted dataset.

For details of the block permutation procedure, see the Supplementary Materials in Ghorbani et al. (2025).

Note that the sim.procedures covers both pure and block permutation methods.

Value

A list of length min(nperm, factorial(nblocks) - 1). Each element is a matrix with the same number of columns as X; the third column contains the block-permuted time values.

References

Ghorbani, M., Vafaei, N. and Myllymäki, M. (2025). A kernel-based test for the first-order separability of spatio-temporal point processes, TEST.

See Also

sim.procedures

Examples

set.seed(123)
X <- cbind(runif(100), runif(100), runif(100, 0, 10))
perms <- block.permut(nblocks = 5, X = X, nperm = 10)
head(perms[[1]], 5)

Compute bandwidths and edge-correction factors for spatio-temporal kernel intensity estimation

Description

Computes spatial and temporal bandwidths for kernel-based estimation of the intensity of a spatio-temporal point pattern. The spatial bandwidth is estimated from the spatial coordinates using Diggle's (1985) mean-square error method via bw.diggle. The temporal bandwidth is estimated from the time coordinates using the Sheather–Jones direct plug-in method (bw = "SJ-dpi") as implemented in density.

Usage

calc.bandwidths.and.edgecorr(
  X,
  s.region,
  t.region,
  n.grid,
  epsilon = NULL,
  delta = NULL
)

Arguments

Details

Edge-correction factors are computed for Gaussian kernels as the kernel mass inside the observation window: c(x) = \int_W K_h(u-x)\,du. The temporal correction is computed exactly over t.region. The spatial correction is computed using the Diggle edge-correction method as implemented in spatstat.explore, which accounts for the polygonal spatial observation window.

The spatial window is converted to an owin object, and the spatial bandwidth is estimated using Diggle's mean-square error method via bw.diggle applied to the corresponding ppp object. Spatial edge-correction factors are computed using the Diggle edge-correction method implemented in spatstat.explore.

The temporal bandwidth is selected using the Sheather–Jones direct plug-in method (bw = "SJ-dpi") as implemented in density over t.region.

The returned edge-correction factors are kernel masses inside the window. If you use them for intensity estimation with edge correction, typical usage is to divide by these factors.

Value

A list with components:

bw

Numeric vector of length 3: c(epsilon, epsilon, delta).

time

Numeric vector of length nrow(X) giving temporal edge masses \int_{tmin}^{tmax} K_\delta(u-t_i)\,du.

space

Numeric vector of length nrow(X) giving spatial edge-correction factors computed using the Diggle method for the polygonal window s.region.

References

Baddeley A, Rubak E, Turner R (2015). Spatial Point Patterns: Methodology and Applications with R. Chapman and Hall/CRC Press.

Diggle, P.J. (1985). A Kernel Method for Smoothing Point Process Data. Journal of the Royal Statistical Society, Series C, 34, 138–147.

See Also

bw.diggle, density

Examples

  set.seed(123)
  X <- cbind(runif(100), runif(100), runif(100, 0, 10))
  s.region <- matrix(c(0,0, 1,0, 1,1, 0,1), ncol = 2, byrow = TRUE)
  t.region <- c(0, 10)
  n.grid <- c(25, 25, 20)
  res <- calc.bandwidths.and.edgecorr(X, s.region, t.region, n.grid)
  str(res)

Validate common arguments for spatio-temporal grid-based routines

Description

Checks the validity of common inputs used by spatio-temporal estimation and simulation routines: event locations X, spatial window specification s.region, temporal window t.region, and grid resolution n.grid.

Usage

check.args(X, s.region, t.region, n.grid = c(25L, 25L, 20L))

Arguments

Details

This function is intended as a lightweight argument checker used internally by multiple functions. It can also be called directly for debugging input issues.

Value

Invisibly returns NULL. Called for its side effect of throwing an error if inputs are invalid.

Examples

X <- cbind(runif(100), runif(100), runif(100, 0, 10))
s.region <- matrix(c(0,0, 1,0, 1,1, 0,1), ncol = 2, byrow = TRUE)
t.region <- c(0, 10)
check.args(X, s.region, t.region, n.grid = c(25, 25, 20))

Chi-squared test for first-order separability of a spatio-temporal point process

Description

Performs a chi-squared test for testing first-order separability of a spatio-temporal point process. Two procedures are available:

"pure_per"

Classical asymptotic chi-squared test of independence on a space–time count table.

"block_per"

Monte Carlo permutation test based on block-wise permutations of the time component.

Usage

chi2.test(
  X,
  sim.procedure = c("pure_per", "block_per"),
  nblocks = 5L,
  nperm = 199L,
  n.time = 2L,
  n.space = 3L,
  t.region = c(0, 1),
  s.region = c(0, 1, 0, 1)
)

Arguments

Details

The classical procedure ("pure_per") applies a chi-squared test of independence to the n.space^2 by n.time contingency table of counts.

The permutation procedure ("block_per") generates up to nperm block-permuted datasets under the null using sim.procedures with method = "block", recomputes the chi-squared statistic for each, and returns a Monte Carlo p-value computed as (1 + \#\{T_i \ge T_{obs}\})/(nperm + 1).

Value

Numeric scalar: the p-value of the test.

Author(s)

Mohammad Ghorbani mohammad.ghorbani@slu.se
Nafiseh Vafaei nafiseh.vafaei@slu.se

References

Ghorbani M., Vafaei N., Dvořák J., Myllymäki M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics and Data Analysis, 161, 107245.

Ghorbani, M., Vafaei, N. and Myllymäki, M. (2025). A kernel-based test for the first-order separability of spatio-temporal point processes, TEST.

See Also

chisq.test.stPP, sim.procedures, block.permut

Examples

set.seed(124)
lambda <- get.lambda.function(N = 200, g = 50, model = 4)
Lmax <- get.lambda.max(N = 200, g = 50, model = 4)
X <- rstpoispp(lambda, Lmax)


# Classical chi-squared test
chi2.test(X, sim.procedure = "pure_per", n.time = 2, n.space = 3)

# Monte Carlo permutation test with blocks
chi2.test(X, sim.procedure = "block_per", nblocks = 5, nperm = 100)

Chi-squared test for first-order separability of a spatio-temporal point process

Description

Performs the classical (asymptotic) chi-squared test of first-order separability by constructing a space-time contingency table of counts and applying a chi-squared test of independence.

Usage

chisq.test.stPP(
  X,
  n.space = 2L,
  n.time = 3L,
  s.region = c(0, 1, 0, 1),
  t.region = c(0, 1)
)

Arguments

Details

The spatial domain is partitioned into n.space bins in each coordinate direction (yielding n.space^2 spatial cells), and the temporal domain is partitioned into n.time intervals. Bin boundaries are defined using empirical quantiles of the observed coordinates, with the first/last boundaries fixed to the provided spatial and temporal windows.

Events falling outside s.region or t.region are ignored (with a warning). If the data contain many ties, quantile-based boundaries may coincide; in that case reduce n.space/n.time or jitter the coordinates slightly.

This implementation uses chisq.test on the contingency table of space–time counts. If expected counts are very small, the chi-squared approximation may be poor; in that case consider using a Monte Carlo approach (e.g., block permutation) as implemented in chi2.test.

Value

A list with components:

chisq_s

Numeric scalar. The chi-squared test statistic.

chisq_p

Numeric scalar. The p-value of the chi-squared test.

counts

Integer matrix of dimension n.space^2 by n.time containing the space–time counts.

Author(s)

Jiří Dvořák dvorak@karlin.mff.cuni.cz

References

Ghorbani M., Vafaei N., Dvořák J., Myllymäki M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics and Data Analysis, 161, 107245.

See Also

chi2.test, chisq.test

Examples

lambda <- get.lambda.function(N = 200, g = 50, model = 4)
Lmax <- get.lambda.max(N = 200, g = 50, model = 4)
X <- rstpoispp(lambda, Lmax)
result <- chisq.test.stPP(X, n.space = 2, n.time = 2)
print(result)

dHSIC test for first-order separability of a spatio-temporal point process

Description

Performs a nonparametric test of first-order separability between space and time in a spatio-temporal point process using the d-variable Hilbert–Schmidt independence criterion (dHSIC). The test statistic evaluates whether the spatio-temporal intensity \rho(u,t) can be written in the multiplicative form \rho_1(u)\rho_2(t),\,\, (u,t)\in \mathbb{R}^2\times\mathbb{R}, where \rho_1 and \rho_2 are nonnegative measurable functions.

Usage

dHS.test(
  X,
  sim.procedure = c("pure_per", "block_per"),
  nblocks = 7L,
  nperm = 1999L,
  nsim = 199L,
  bandwidth = NULL
)

Arguments

Details

Two permutation strategies are supported:

"pure_per"

Randomly permutes the time component across events.

"block_per"

Uses block-wise permutation of the time component via sim.procedures with method = "block" to preserve short-range temporal dependence.

The Monte Carlo p-value is computed with the standard +1 correction: (1 + \#\{T_i \ge T_{obs}\})/(B + 1), where B is the number of permutations.

Value

A list with components:

p.value

Monte Carlo p-value based on the adaptive-bandwidth Gaussian kernel.

p.value.bw

Monte Carlo p-value based on the fixed-bandwidth Gaussian kernel, or NA if bandwidth = NULL.

bandwidth_data

Bandwidth selected by dHSIC::dhsic(..., kernel = "gaussian").

Note

The dHSIC method is implemented via the dHSIC package. When sim.procedure = "pure_per", dHS.test() internally calls dhsic from R package dHSIC for computational efficiency.

References

Ghorbani, M., Vafaei, N. and Myllymäki, M. (2025). A kernel-based test for the first-order separability of spatio-temporal point processes, TEST.

See Also

sim.procedures, block.permut, chi2.test

Examples

if (requireNamespace("dHSIC", quietly = TRUE)) {
  set.seed(123)
  X <- cbind(runif(100), runif(100), runif(100, 0, 10))

  # Pure permutation test
  result <- dHS.test(sim.procedure = "pure_per",
                     X = X, nsim = 199, bandwidth = 0.05)
  print(result$p.value)

  # Block permutation test
  result_block <- dHS.test(sim.procedure = "block_per", X = X,
                           nblocks = 5, nperm = 100, bandwidth = 0.05)
  print(result_block$p.value.bw)
}

Kernel-based intensity estimation on a space-time grid and its components, and test statistics for first-order separability

Description

Computes kernel-smoothed estimates of spatial, temporal, separable, and non-separable spatio-temporal intensity functions on a regular space-time grid, together with separability diagnostics used in first-order separability analysis.

Usage

estimate.intensity.pixel(X, s.region, t.region, n.grid, edge, owin = NULL)

Arguments

Details

The estimator uses product Gaussian kernels with supplied bandwidths and (Gaussian) edge-correction factors, typically produced by calc.bandwidths.and.edgecorr.

Value

A list with grid coordinates x,y,t, intensity estimates, the diagnostic S.fun, its marginal summaries S.space and S.time, and deviation measures.

References

Ghorbani, M., Vafaei, N., Dvořák, J., and Myllymäki, M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics & Data Analysis, 161, 107245.

See Also

S.based.functions, calc.bandwidths.and.edgecorr

Examples

n <- 100
X <- cbind(x = stats::runif(n), y = stats::runif(n), t = stats::runif(n, 0, 10))
s.region <- matrix(c(0,0, 1,0, 1,1, 0,1), ncol=2, byrow=TRUE)
t.region <- c(0, 10)
n.grid <- c(10, 10, 5)
edge <- list(bw = c(0.05, 0.05, 0.5), space = rep(1, n), time = rep(1, n))
res <- estimate.intensity.pixel(X, s.region, t.region, n.grid, edge)
str(res)

Kernel intensity estimates of a spatio-temporal point process at observed points and its components, and test statistics for first-order separability

Description

Computes kernel-based spatial, temporal, separable, and non-separable intensity estimates evaluated at the observed spatio-temporal event locations. The function also returns the separability diagnostic S_i and global deviation measures quantifying departures from first-order separability.

Usage

estimate.intensity.point(X, n.grid, edge)

Arguments

Details

Pairwise Gaussian kernel weights are computed in each dimension and diagonal entries are set to zero to remove self-contributions.

Value

A list with components S.fun, deviation measures, and estimated intensity components at the observed points.

See Also

dnorm

Examples

X <- cbind(stats::runif(50), stats::runif(50), stats::runif(50))
edge <- list(bw = c(0.1, 0.1, 0.1), space = 1, time = 1)
res <- estimate.intensity.point(X, n.grid = 50, edge = edge)
str(res)

Kernel Estimation of the Spatio-Temporal Intensity Function and Its Components, and Test Statistics for First-Order Separability

Description

It returns estimates of the full spatio-temporal intensity together with its spatial and temporal marginal components. The output also includes the test statistic S(u,t) and its spatial and temporal marginal profiles S_{\mathrm{space}}(u) and S_{\mathrm{time}}(t) (equations (11)–(13) in Ghorbani et al., 2021), as well as several deviation tests (e.g., equation (15) in Ghorbani et al., 2021) that quantify departures from first-order separability.

Usage

estimate.st.intensity(
  X,
  s.region,
  t.region,
  at = c("pixels", "points"),
  n.grid = c(25, 25, 20)
)

Arguments

Details

The estimation follows the kernel framework of Ghorbani et al. (2021). A Gaussian kernel is applied in each spatial and temporal dimension. Spatial bandwidths are estimated using Diggle’s method, and the temporal bandwidth is obtained by the Sheather–Jones direct plug-in (SJ-DPI) approach (see calc.bandwidths.and.edgecorr for implementation details). Edge corrections and bandwidths are computed using the calc.bandwidths.and.edgecorr.

The non-separable intensity estimate is

\hat{\rho}(u,t) = \sum_{i=1}^n \frac{K_\epsilon^2\!\left(u-u_i\right)}{C_{W,\epsilon}(u_i)} \frac{K_\delta^1\!\left(t-t_i\right)}{C_{T,\delta}(t_i)},

where k_b(v) = k(v/b)/b^d is a d-dimensional kernel with bandwidth b > 0, and K^1 and K^2 are Gaussian kernels with spatial and temporal bandwidths \epsilon and \delta. C_{W,\epsilon}(u_i) and C_{T,\delta}(t_i) are spatial and temporal edge corrections, respectively. For estimate of the separable counterparts and more details see equations (3)-(5) in Ghorbani et al., 2021).

The separability test is:

S(u, t) = \frac{\hat \rho(u, t)}{\hat\rho_{\mathrm{space}}(u)\hat\rho_{\mathrm{time}}(t)/n},

where n is the number of observed points.

Several deviation statistics are provided to quantify departures from separability

• deviation.t1: Integral of absolute deviations \int_{W\times T}|\hat\rho(u,t) - \hat\rho_{\mathrm{sep}}(u,t)|\mathrm{d}u\mathrm{d}t.

• deviation.t2: Integral of absolute deviation of inverse intensities \int_{W\times T}|1/\hat\rho(u,t) - 1/\hat\rho_{\mathrm{sep}}(u,t)|\mathrm{d}u\mathrm{d}t.

• deviation.t3: Sum of log-ratio deviations \sum_i\big(\log(\hat\rho(u_i,t_i)- \log(\hat\rho_{\mathrm{sep}}(u_i,t_i))\big).

• deviation.t4: Integral of the S-function.

When at = "points", intensities and diagnostics are evaluated at the observed event locations.

Value

A list containing:

x, y, t

Grid vectors for x, y, and t (returned only when at = "pixels").

epsilon

Estimated spatial bandwidth used in Gaussian kernels.

delta

Estimated temporal bandwidth.

SPat.intens

Estimated spatial intensity surface (pixels only).

TeM.intens

Estimated temporal intensity profile (pixels only).

sep.intens

Separable spatio-temporal intensity estimate.

nonsep.intens

Non-separable spatio-temporal intensity estimate.

S.fun

Test function S(u,t) as the ratio of non-separable to separable intensity estimates.

S.space

Marginal sum of S(u,t) over time (space profile).

S.time

Marginal sum of S(u,t) over space (time profile).

deviation.t1

Deviation statistic: Integral of absolute deviations.

deviation.t2

Deviation statistic: absolute deviation of inverse intensities.

deviation.t3

Deviation statistic: sum of log-ratio deviations.

deviation.t4

Total integral of the S-function.

Author(s)

Mohammad Ghorbani mohammad.ghorbani@slu.se
Nafiseh Vafaei nafiseh.vafaei@ltu.se

References

Ghorbani M., Vafaei N., Dvořák J., Myllymäki M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics & Data Analysis, 161, 107245.

See Also

calc.bandwidths.and.edgecorr, global.envelope.test, chi2.test, dHS.test

Examples

set.seed(123)
X <- cbind(runif(100), runif(100), runif(100, 0, 10))
s.region <- matrix(c(0, 0, 1, 0, 1, 1, 0, 1), ncol = 2, byrow = TRUE)
t.region <- c(0, 10)
result <- estimate.st.intensity(X, s.region, t.region, at = "pixels", n.grid = c(64, 64, 32))
str(result$S.fun)
image(result$S.fun[,,5], main = "S-function slice at t[5]",
       col = topo.colors(50))
contour(result$S.fun[,,5], add=TRUE)
str(result[c("deviation.t1","deviation.t3","deviation.t4")])

Construct spatio-temporal intensity functions with controlled separability

Description

Returns an intensity function \lambda(x,y,t) corresponding to one of four models used for simulation experiments in Ghorbani et al. (2021). Each model is a mixture of a separable "background" component and a structured (generally non-separable) spatio-temporal Gaussian bump. The models provide different degrees of space–time separability, allowing for controlled experiments on separability testing.

Usage

get.lambda.function(
  N,
  g,
  model = 1L,
  mu1 = 0.5,
  sd1 = 0.2,
  mu2 = c(0.5, 0.5),
  sd2 = c(0.2, 0.2),
  mu3 = c(0.3, 0.3, 0.2),
  sd3 = c(0.05, 0.05, 0.05)
)

Arguments

Details

The returned function is intended for use in simulation (e.g., for generating spatio-temporal Poisson point patterns under varying degrees of separability).

The intensity is constructed as:

\lambda(x,y,t) = \lambda_{\mathrm{bg}}(x,y,t) + g\, f_{st}(x,y,t),

where f_{st} is a nonnegative 3D Gaussian density (via norm3d) and the background term \lambda_{\mathrm{bg}} depends on model:

1

Homogeneous background: \lambda_{\mathrm{bg}}(x,y,t) = (N-g)

2

Temporal inhomogeneity only: \lambda_{\mathrm{bg}}(x,y,t) = (N-g)\, f_t(t), where f_t is a 1D Gaussian density (dnorm).

3

Spatial inhomogeneity only: \lambda_{\mathrm{bg}}(x,y,t) = (N-g)\, f_s(x,y), where f_s is a 2D Gaussian density (norm2d).

4

Separable spatial-temporal inhomogeneity: \lambda_{\mathrm{bg}}(x,y,t) = (N-g)\, f_s(x,y)\, f_t(t).

Note: since Gaussian densities can exceed 1 for small standard deviations, N is best interpreted as a scaling parameter used by the calling simulator. Ensure \lambda(x,y,t) is nonnegative over the intended domain.

See more details in Ghorbani et al. (2021), Section 6.1.

Value

A function of the form function(x, y, t) representing the selected intensity surface.

Note

This function is primarily intended for generating intensity functions used in simulation studies. In particular, rstpoispp calls get.lambda.function() internally to construct intensity models for simulating spatio-temporal Poisson point processes with controlled separability.

Author(s)

Mohammad Ghorbani mohammad.ghorbani@slu.se
Nafiseh Vafaei nafiseh.vafaei@slu.se

References

Ghorbani, M., Vafaei, N., Dvořák, J., and Myllymäki, M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics & Data Analysis, 161, 107245.

See Also

rstpoispp, norm3d, norm2d, chi2.test, dHS.test

Examples

# Choose model 4: non-separable spatio-temporal intensity
lambda <- get.lambda.function(N = 210, g = 50, model = 4)
lambda(0.5, 0.5, 0.5)  # Evaluate intensity at center of space-time domain

# Visualize spatial intensity at fixed time for model 2
lambda2 <- get.lambda.function(N = 200, g = 50, model = 2)
x <- y <- seq(0, 1, length.out = 100)
z <- outer(x, y, function(x, y) lambda2(x, y, t = 0.5))

fields::image.plot(
x, y, z,
main = "Intensity at t = 0.5 (Model 2)",
col = topo.colors(50))

Upper bound for spatio-temporal intensity models

Description

Computes a practical upper bound for the spatio-temporal intensity models used in get.lambda.function. The bound is intended for thinning/rejection sampling in simulation routines such as rstpoispp.

Usage

get.lambda.max(
  N,
  g,
  model = 1L,
  mu1 = 0.5,
  sd1 = 0.2,
  mu2 = c(0.5, 0.5),
  sd2 = c(0.2, 0.2),
  mu3 = c(0.3, 0.3, 0.2),
  sd3 = c(0.05, 0.05, 0.05)
)

Arguments

Details

The bound is computed using analytic maxima of Gaussian density components (at their modes), which yields a conservative and fast-to-evaluate upper bound when the component functions are Gaussian product densities.

The intensity models are mixtures of a background term and a structured spatio-temporal Gaussian bump. This function returns an upper bound obtained by evaluating each Gaussian density component at its mode. This upper bound is typically sufficient for rejection sampling when generating realizations of inhomogeneous Poisson or Cox point processes.

If norm2d and norm3d in this package are Gaussian product densities (independent components), then the maxima are available in closed form:

• \max_t \phi(t;\mu,\sigma) = 1/(\sigma\sqrt{2\pi})

• \max_{x,y} \phi(x;\mu_x,\sigma_x)\phi(y;\mu_y,\sigma_y) = 1/(2\pi\sigma_x\sigma_y)

• \max_{x,y,t} \prod_{d=1}^3 \phi(\cdot;\mu_d,\sigma_d) = 1/((2\pi)^{3/2}\sigma_x\sigma_y\sigma_t)

If your norm2d/norm3d use a different parameterization, this bound should be updated accordingly.

Value

Numeric scalar. A conservative upper bound for the selected intensity model.

Author(s)

Nafiseh Vafaei nafiseh.vafaei@slu.se
Mohammad Ghorbani mohammad.ghorbani@slu.se

References

Ghorbani, M., Vafaei, N., Dvořák, J., and Myllymäki, M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics and Data Analysis, 161, 107245.

See Also

get.lambda.function, rstpoispp

Examples

# Example 1: Homogeneous model (Model 1)
get.lambda.max(N = 200, g = 50, model = 1)

# Example 2: Non-separable spatio-temporal model (Model 4)
get.lambda.max(N = 200, g = 50, model = 4)

Global envelope test for spatio-temporal separability using S-function

Description

Performs a global envelope test of the null hypothesis of first-order separability for a spatio-temporal point process. The observed separability diagnostics S(u,t), S_{\mathrm{space}}(u) and/or S_{\mathrm{time}}(t) are compared to a reference distribution obtained from permuted versions of the data.

Usage

global.envelope.test(
  X,
  sim.procedure = c("pure_per", "block_per"),
  nsim = 25L,
  nblocks = 5L,
  nperm = 199L,
  n.grid = c(20L, 20L, 10L),
  s.region = matrix(c(0, 0, 1, 0, 1, 1, 0, 1), ncol = 2, byrow = TRUE),
  t.region = c(0, 1),
  owin = NULL,
  eps = NULL,
  del = NULL,
  tests = c("S.test", "S.space.test", "S.time.test"),
  GET.args = NULL
)

Arguments

Details

Two permutation strategies are supported:

"pure_per"

Pure permutation: randomly permutes the time coordinates.

"block_per"

Block permutation: permutes time in blocks to preserve short-range temporal dependence.

The GET package is used to construct global envelopes and compute p-values.

The null hypothesis is

H_0:\ \rho(u,t)=\rho_{\mathrm{space}}(u)\rho_{\mathrm{time}}(t).

The function computes the chosen diagnostics using S.based.functions on a pixel grid (at="pixels") and applies global_envelope_test.

To keep curve lengths identical (required by GET), any NA values induced by an owin mask are removed using the same indices for the observed and all simulated curves.

Value

A list with components:

Bandwidth_s

Spatial bandwidth used.

Bandwidth_t

Temporal bandwidth used.

S.test, S.space.test, S.time.test

For each requested test: p-values for ERL and AREA envelopes, plus optional plots if created.

Author(s)

Nafiseh Vafaei nafiseh.vafaei@slu.se
Mohammad Ghorbani mohammad.ghorbani@slu.se

References

Ghorbani, M., Vafaei, N., Dvořák, J., and Myllymäki, M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics & Data Analysis, 161, 107245.

See Also

S.based.functions, sim.procedures, block.permut, global_envelope_test, plot.global_envelope

Examples

if (requireNamespace("GET", quietly = TRUE)) {
  set.seed(123)
  X <- cbind(stats::runif(100), stats::runif(100), stats::runif(100, 0, 1))
  s.region <- matrix(c(0,0, 1,0, 1,1, 0,1), ncol = 2, byrow = TRUE)
  t.region <- c(0, 1)

  res <- global.envelope.test(
    X = X,
    sim.procedure = "pure_per",
    nsim = 19,
    n.grid = c(10,10,10),
    s.region = s.region,
    t.region = t.region,
    tests = c("S.test","S.time.test")
  )
  str(res)
}

Bivariate normal density with independent components

Description

Evaluates the density of a bivariate normal distribution with mean vector \mu = (\mu_x,\mu_y) and diagonal covariance matrix (independent components). The density is the product of the two univariate normal densities:

f(x,y) = \phi(x;\mu_x,\sigma_x)\,\phi(y;\mu_y,\sigma_y).

Usage

norm2d(x, y, mu = c(0, 0), sd = c(1, 1), log = FALSE)

Arguments

Value

Numeric vector of densities (or log-densities) with length determined by standard recycling rules for x and y.

Author(s)

Mohammad Ghorbani mohammad.ghorbani@slu.se
Nafiseh Vafaei nafiseh.vafaei@slu.se

Examples

# Evaluate the density at the peak
norm2d(0.5, 0.5, mu = c(0.5, 0.5), sd = c(0.2, 0.2))

# Evaluate at multiple x values
norm2d(c(0.3, 0.7), 0.5, mu = c(0.5, 0.5), sd = c(0.2, 0.2))

# Visualize on a grid
x <- y <- seq(0, 1, length.out = 100)
f <- Vectorize(function(x, y) norm2d(x, y, mu = c(0.5, 0.5), sd = c(0.2, 0.2)))
z <- outer(x, y, f)
image(x, y, z, col = terrain.colors(50), main = "Bivariate Normal Intensity")
contour(x, y, z, add = TRUE)

Trivariate normal density with independent components in space-time

Description

Evaluates a trivariate normal density on (x,y,t) with independent components (diagonal covariance). The density is the product of three univariate normal densities:

f(x,y,t) = \phi(x;\mu_x,\sigma_x)\,\phi(y;\mu_y,\sigma_y)\,\phi(t;\mu_t,\sigma_t).

Usage

norm3d(x, y, t, mu = c(0.3, 0.3, 0.2), sd = c(0.05, 0.05, 0.05), log = FALSE)

Arguments

Value

Numeric vector of densities (or log-densities) with length determined by standard recycling rules for x, y, and t.

Author(s)

Mohammad Ghorbani mohammad.ghorbani@slu.se

References

Ghorbani, M., Vafaei, N., Dvořák, J., and Myllymäki, M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics & Data Analysis, 161, 107245.

See Also

norm2d, get.lambda.function, estimate.st.intensity

Examples

norm3d(0.3, 0.3, 0.2)  # peak value at the mean (with default parameters)
norm3d(c(0.2, 0.3), 0.3, 0.2)

x <- y <- seq(0, 1, length.out = 100)
z <- outer(x, y, function(x, y) norm3d(x, y, t = 0.2))
image(x, y, z, col = heat.colors(50), main = "Spatial slice of norm3d at t = 0.2")

Plot original vs permuted time ordering

Description

Produces a side-by-side plot comparing the temporal component of an original spatio-temporal point pattern with that of a permuted (or block-permuted) version. This graphical diagnostic is intended to assess the effect of temporal reordering procedures used in separability tests.

Usage

plot_procedures(
  original,
  permuted,
  title = "Permutation",
  col = c("blue", "red"),
  pch = 1,
  ...
)

Arguments

Details

The function is commonly employed to visualize temporal permutations generated by procedures such as sim.procedures or block.permut, which underpin pure and block permutations-based inference for first-order separability (see Ghorbani et al., 2021, Section 3.2).

The function uses base R graphics to display two panels:

1. The temporal ordering of the original point pattern.

2. The temporal ordering after permutation or block permutation.

This diagnostic is particularly useful when validating permutation-based inference procedures such as chi2.test, dHS.test, and global.envelope.test.

Value

The function is invoked for its side effect of producing a plot and returns no value.

Author(s)

Nafiseh Vafaei nafiseh.vafaei@slu.se

References

Ghorbani, M., Vafaei, N., Dvořák, J., & Myllymäki, M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics & Data Analysis, 161, 107245.

See Also

sim.procedures, block.permut, dHS.test, chi2.test

Examples

set.seed(123)
X <- cbind(x = runif(20), y = runif(20),
           time = seq(0, 1, length.out = 20))

# Example: visualize pure permutation
sim_pure <- sim.procedures(X, nperm = 1, method = "pure")[[1]]
plot_procedures(X, sim_pure, title = "Pure Permutation")

# Example: visualize block permutation (requires Block_permutation)
if (requireNamespace("combinat", quietly = TRUE)) {
  sim_block <- block.permut(nblocks = 4, X = X, nperm = 1)[[1]]
  plot_procedures(X, sim_block, title = "Block Permutation")
}

Plot spatio-temporal determinantal point process (DPP) realizations

Description

Produces diagnostic plots for spatio-temporal determinantal point process (DPP) simulations or fitted models. The function formats the plot title based on the spatial and temporal interaction parameters \alpha_s and \alpha_t, automatically displaying either the parameters themselves or their reciprocals when greater than 1.

Usage

plot_stDPP(data, type = c("3D", "space", "time"), alpha_s, alpha_t)

Arguments

Details

If \alpha_s > 1 and \alpha_t > 1, the title displays \alpha_s^{-1} and \alpha_t^{-1}, which correspond to interaction ranges. Otherwise, the parameters are shown directly.

The function then calls plot.ST.pp() to produce the actual plot.

Value

No return value. The function is called for its side effect of producing a diagnostic plot.

Examples

# Simulate a stationary separable Matérn ST-DPP
sim <- rstDPP(
  mode     = "stationary",
  model    = "S",
  spectral = "matern",
  alpha_s  = 10,
  alpha_t  = 4.7,
  nu       = 2,
  eps      = 1,
  lambda_max = 70,
  grid_size  = 1.5
)
 plot_stDPP(sim, type = "3D",  alpha_s = 10, alpha_t = 4.7)

Plot spatio-temporal log-Gaussian Cox process (LGCP) realizations as time-sliced maps

Description

Produces a sequence of spatial raster maps displaying the evolution of a simulated or fitted spatio-temporal log-Gaussian Cox process (LGCP) over time. Each map shows the latent Gaussian random field (intensity surface) at a given time slice, with the corresponding observed point events overlaid. This visualization helps interpret temporal evolution and spatial clustering patterns in simulated or fitted LGCP models.

Usage

plot_stlgcp(data)

Arguments

Details

The function plots up to 10 evenly spaced time slices from the latent intensity field and overlays the corresponding point events accumulated up to each time point. Each panel represents a spatial realization at a fixed time t_k, providing a visual summary of the dynamic spatio-temporal structure of the LGCP.

This approach follows the visualization principles used in Ghorbani et al. (2021, 2025), where LGCPs are employed to assess the behavior of separability diagnostics and kernel-based tests in complex, non-separable point process settings. The raster intensity surfaces illustrate latent heterogeneity, while the overlaid points display observed event clustering relative to the underlying field.

Time slices are selected at evenly spaced quantiles of the temporal domain, and each plot includes:

• A raster map of the latent intensity surface at time t_k;

• White contour lines showing equal-intensity regions;

• Overlaid black points representing events observed up to t_k.

The resulting plots are arranged into a single grid layout using the patchwork package for ease of comparison.

Value

A combined ggplot object displaying up to ten spatial raster maps arranged in a grid layout (by default, two rows and up to five columns). Each panel corresponds to one time slice. The function returns the combined plot object.

Author(s)

Nafiseh Vafaei nafiseh.vafaei@slu.se

References

Ghorbani M., Vafaei N., Dvořák J., Myllymäki M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics & Data Analysis, 161, 107245.

See Also

Gauss.st.F for simulating spatio-temporal Gaussian random fields; get.lambda.function and rstpoispp for generating intensity-based spatio-temporal point processes.

Examples

# Example: visualize a spatio-temporal LGCP simulation
out <- rstLGCPP(xlim = c(0,1),
                ylim = c(0,1),
                tlim = c(0,1),
                grid = c(15,15,10))
plot_stlgcp(data = out)

Plot a spatio-temporal point pattern

Description

Provides flexible visualization tools for spatio-temporal point patterns. Depending on the chosen display type, the function produces one of three plots:

• A 3D scatterplot showing event locations in space–time ("3D");

• A 2D spatial projection of points in the spatial plane ("space");

• A temporal histogram with an overlaid kernel density curve ("time").

The input dataset must contain three columns corresponding to the spatial (x, y) and temporal (t) coordinates of events.

Usage

plot_stpp(data, type = c("3D", "space", "time"), time_bins = 30, title = NULL)

Arguments

Details

The function serves as an exploratory tool for investigating the spatial, temporal, or joint space–time structure of point pattern data. Such visualization is often a first step before conducting statistical analyses of separability or intensity modeling (see Ghorbani et al., 2021).

3D mode

Displays events in a 3D coordinate system using the scatterplot3d package, allowing a quick assessment of clustering and temporal trends in space–time.

Spatial mode

Projects points onto the spatial plane (x–y), showing spatial structure independent of time.

Temporal mode

Displays the marginal temporal distribution of events as a histogram with a kernel density overlay, facilitating the visual detection of temporal nonstationarity.

Visualization is particularly useful for validating the realism of simulated data (e.g. from rstpoispp or Gauss.st.F) and for preliminary inspection prior to applying tests such as chi2.test, global.envelope.test, or dHS.test.

Value

Produces a plot as a side effect. Nothing is returned.

Note

The 3D visualization requires the scatterplot3d package.

Author(s)

Nafiseh Vafaei nafiseh.vafaei@slu.se

References

Ghorbani M., Vafaei N., Dvořák J., Myllymäki M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics & Data Analysis, 161, 107245.

See Also

rstpoispp, estimate.st.intensity, plot_stlgcp, chi2.test, dHS.test

Examples

set.seed(123)
X <- cbind(runif(100), runif(100), runif(100, 0, 10))

# Visualize point pattern in 3D space–time
plot_stpp(X, type = "3D")

# View spatial projection
plot_stpp(X, type = "space")

# Inspect temporal distribution
plot_stpp(X, type = "time", time_bins = 20)

Apply a circular random shift to the temporal component of a spatio-temporal point pattern

Description

Performs a circular random shift of the temporal coordinate in a spatio-temporal point pattern. This operation preserves the spatial configuration while randomizing the temporal component under the assumption of temporal stationarity. The shift amount is drawn uniformly from [0, 1] and applied modulo 1, ensuring that the time window is maintained. The resulting dataset can be used to construct null models for hypothesis testing of first-order separability or temporal independence.

Usage

random.shift(X, shifted_col = 3)

Arguments

Details

The circular random shift is a common resampling procedure for generating null models of temporal randomness while preserving the overall temporal marginal distribution and spatial structure.

For each dataset, a single uniform random shift value \Delta \sim \mathrm{Uniform}(0,1) is drawn and added to the temporal coordinate. The shifted times are then wrapped around the unit interval:

t_i^{*} = (t_i + \Delta) \bmod 1, \quad i = 1, \dots, n.

Value

A matrix or data frame (matching the input type) with the time column shifted modulo 1.

References

Ghorbani, M., Vafaei, N., Dvořák, J., and Myllymäki, M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics & Data Analysis, 161, 107245.

Examples

set.seed(123)
X <- cbind(runif(100), runif(100), runif(100))  # x, y, t
X_shifted <- random.shift(X)

# Compare original and shifted time values
head(X[,3])
head(X_shifted[,3])

# Verify shift visually
plot(X[, 3], type = "o", col = "blue", ylab = "Time", xlab = "Index",
     main = "Original vs Shifted Times")
lines(X_shifted[, 3], type = "o", col = "red")
legend("topright", legend = c("Original", "Shifted"),
       col = c("blue", "red"), lty = 1, pch = 1)

Simulate a spatio-temporal Determinantal Point Process (DPP) based on spectral density

Description

Generates a realization of a spatio-temporal determinantal point process (DPP) using a user-defined spectral density model. The function supports both separable (model = "S") and non-separable (model = "NS") dependence structures and allows for either exponential or Mat\'ern-type spectral densities. The simulation is performed on a 3D spatio-temporal frequency grid and can incorporate user-specified intensity functions for thinning.

Usage

rstDPP(
  mode = c("stationary", "inhom"),
  model = c("S", "NS"),
  spectral = c("exp", "matern"),
  alpha_s,
  alpha_t = NULL,
  lambda_max,
  nu = 2,
  eps = 1,
  lambda_s = NULL,
  lambda_non_s = NULL,
  grid_size = 4
)

Arguments

Details

This function implements a spectral simulation method for spatio-temporal DPPs, following the theoretical framework introduced in Vafaei et al. (2023) and Ghorbani et al. (2025).

The algorithm proceeds as follows:

1. Construct a 3D grid of spatial and temporal frequency components (\omega_x, \omega_y, \tau).

2. Evaluate the chosen spectral density \phi(\omega, \tau) across the grid.

3. Use the resulting spectral values as eigenvalues to simulate a realization of a DPP via spatstat.model::rdpp().

4. Optionally apply thinning using a user-defined intensity function \lambda(u, t), scaled by lambda_max, to induce inhomogeneity.

Two spectral families are supported:

• Exponential form:

  \phi(\omega, \tau) \propto \exp\left[-(\pi \alpha_s |\omega|)^2\right] \left(1 + 4 (\pi \alpha_t \tau)^2\right)^{-1}.

• Mat\'ern-type form:

  \phi_{\epsilon}(\omega, \tau) \propto \left( \alpha_s^2 \alpha_t^2 + \alpha_t^2 |\omega|^2 + \alpha_s^2 \tau^2 + \epsilon |\omega|^2 \tau^2 \right)^{-\nu},

  where \epsilon \in [0, 1] determines the degree of separability between space and time.

This framework enables simulation of spatio-temporal point patterns that exhibit varying degrees of spatial–temporal dependence, providing a versatile tool for evaluating separability tests and modeling non-separable dynamics.

Value

A numeric matrix with three columns (x, y, t) representing the retained spatio-temporal events after thinning.

References

Vafaei, N., Ghorbani, M., Ganji, M., and Myllymäki, M. (2023). Spatio-temporal determinantal point processes. arXiv:2301.02353.

Ghorbani, M., Vafaei, N., and Myllymäki, M. (2025). A kernel-based test for the first-order separability of spatio-temporal point processes. TEST, 34, 580-611. https://doi.org/10.1007/s11749-025-00972-y

See Also

plot_stpp for visualizing spatio-temporal point patterns.

Examples

# Simulate a stationary separable Mat\'ern ST-DPP
 if (requireNamespace("spatstat", quietly = TRUE)) {
sim <- rstDPP(
  mode     = "stationary",
  model    = "S",
  spectral = "matern",
  alpha_s  = 10,
  alpha_t  = 4.7,
  nu       = 2,
  eps      = 1,
  lambda_max = 70,
  grid_size  = 2
)
 plot_stDPP(sim, type = "3D",  alpha_s = 10, alpha_t = 4.7)
# example 2
# Generate realization
sim <- rstDPP(mode = "stationary",
                model = "S",
                spectral = "matern",
                alpha_s = 10, alpha_t = 4.7,
                nu = 2,
               eps = 1,
               lambda_s = 70,
               lambda_non_s = NULL,
               grid_size = 2,
               lambda_max=70)
head(sim)
}

Simulate a spatio-temporal Log-Gaussian Cox process (LGCP)

Description

Generates a realization of a spatio-temporal LGCP over a user-defined domain. The process is simulated using a log-Gaussian random field combined with a deterministic trend function, and points are generated by thinning a homogeneous Poisson process.

Usage

rstLGCPP(
  xlim = NULL,
  ylim = NULL,
  tlim = NULL,
  grid = c(15, 15, 10),
  mu = NULL,
  Lambda = NULL,
  Lmax = NULL,
  par1 = c(1, 0.05),
  par2 = c(1, 0.06),
  sigmas = c(0.5, 0.5, 1),
  mu_par = c(1.2, 0.25, 5)
)

Arguments

Value

A list with:

st.lgcp

A data frame of simulated spatio-temporal points.

RF

The latent Gaussian field output from Gauss.st.F.

Examples

out <- rstLGCPP(xlim = c(0,1),
                ylim = c(0,1),
                tlim = c(0,1),
                grid = c(15,15,10))
plot_stlgcp(data = out)
plot_stpp(data = out$st.lgcp, type = "3D")

Simulate an inhomogeneous spatio-temporal Poisson point process

Description

Generates a realization of an inhomogeneous Poisson point process (STPP) in space and time using the standard thinning method. The user provides an intensity function \lambda(u, t) and an upper bound L_{\max} on its value over the observation window. The algorithm first samples candidate events uniformly over space and time and then retains each candidate with probability proportional to its normalized intensity \lambda(u,t)/L_{\max}.

Usage

rstpoispp(
  lambda,
  Lmax,
  s.region = splancs::as.points(c(0, 1, 1, 0), c(0, 0, 1, 1)),
  t.region = c(0, 1)
)

Arguments

Details

The method implements the classical thinning algorithm for simulating inhomogeneous Poisson processes:

1. Draw N^* \sim \mathrm{Poisson}(L_{\max} \, |W| \, |T|), where |W| and |T| denote the spatial and temporal window measures.

2. Generate N^* candidate points uniformly over W \times T.

3. Retain each point (u_i, t_i) independently with probability p_i = \lambda(u_i, t_i) / L_{\max}.

The result is a realization of an inhomogeneous STPP with intensity function \lambda(u, t).

This simulator underpins the spatio-temporal framework introduced in Ghorbani et al. (2021, 2025) for studying first-order separability. By selecting appropriate intensity functions (see get.lambda.function), users can generate fully separable, partially separable, or non-separable spatio-temporal patterns, enabling direct evaluation of separability tests such as chi2.test, global.envelope.test, or dHS.test.

Value

A numeric matrix with three columns (x, y, t) representing the retained points from the inhomogeneous Poisson process.

Note

The intensity function \lambda(u, t) should return non-negative numeric values and be bounded above by Lmax across the observation domain.

Author(s)

Mohammad Ghorbani mohammad.ghorbani@slu.se
Nafiseh Vafaei nafiseh.vafaei@ltu.se

References

Ghorbani M., Vafaei N., Dvořák J., Myllymäki M. (2021). Testing the first-order separability hypothesis for spatio-temporal point patterns. Computational Statistics & Data Analysis, 161, 107245.

Ghorbani, M., Vafaei, N. and Myllymäki, M. (2025). A kernel-based test for the first-order separability of spatio-temporal point processes, TEST .

See Also

get.lambda.function to construct spatio-temporal intensity models; get.lambda.max to compute intensity maxima; estimate.st.intensity for intensity estimation; plot_stpp for visualization.

Examples

# Example 1: Simulate a separable spatio-temporal Poisson process
lambda <- get.lambda.function(N = 200, g = 50, model = 1)
Lmax   <- get.lambda.max(N = 200, g = 50, model = 1)
X <- rstpoispp(lambda, Lmax)
head(X)

# Example 2: Non-separable model (Model 4)
lambda <- get.lambda.function(N = 200, g = 50, model = 4)
Lmax   <- get.lambda.max(N = 200, g = 50, model = 4)
sim_data <- rstpoispp(lambda, Lmax)

# Spatial projection of simulated events
plot(sim_data[, 1:2], asp = 1, main = "Spatial Projection of Simulated stPP")
# Example 3: 3D visualization using plot_ST_pp()
plot_stpp(X, type = "3D", title="Realisation of a stPP")

Spatio-temporal analysis of the 2001 FMD outbreak (illustrative example)

Description

This function tests first-order separability for the foot-and-mouth disease (FMD) outbreak data observed at the North Cumbria study region provided by the stpp package.

Usage

run_fmd_example()

Details

The example includes kernel-based estimates of the spatio-temporal intensity and related separability diagnostics on a regular spatio-temporal grid for fmd data as well as the dHSIC test.

The datasets fmd and northcumbria are loaded from the stpp package and are used purely for illustration. They do not contain confidential information and have no scientific value.

Value

An object of class "test" returned by dHS.test.

References

Diggle, P., Rowlingson, B. and Su, T. (2005). Point process methodology for on-line spatio-temporal disease surveillance. Environmetrics, 16, 423–434.

See Also

fmd, northcumbria, S.based.functions, dHS.test

Examples

if (requireNamespace("stpp", quietly = TRUE) ) {

  library(stpp)

  # Load data
  data("fmd", "northcumbria", package = "stpp")

  # Define spatial and temporal region
  s.region <- northcumbria / 1000
  t.region <- c(0, 200)

  # Create a spatio-temporal point pattern
  X <- as.3dpoints(
    fmd[, 1] / 1000,
    fmd[, 2] / 1000,
    fmd[, 3]
  )

  # Define observation window
  ObsW <- spatstat.geom::owin(
    poly = list(
      x = s.region[, 1],
      y = s.region[, 2]
    )
  )

  # Run Hilbert-Schmidt independence test
  oldpar <- par(no.readonly = TRUE)  # save current graphics settings
  on.exit(par(oldpar), add = TRUE)   # restore after example

  dHS.test(
    X,
    sim.procedure = "block_per",
    nblocks = 5L,
    nperm = 1999L,
    nsim = 199L,
    bandwidth = NULL
  )
 }

Generate permuted versions (pure or block) of a spatio-temporal point pattern

Description

Implements two types of permutation procedures for resampling the time component of spatio-temporal point process data:

"pure"

Pure random permutation of the time coordinates.

"block"

Block permutation where the time dimension is divided into consecutive blocks, and permutations are applied at the block level.

These procedures are used for generating surrogate datasets under the null hypothesis of first-order separability.

Usage

sim.procedures(X, nperm = 1999, nblocks = 4, method = c("block", "pure"))

Arguments

Value

A list of nperm matrices. Each matrix is a permuted version of the original input X, where the third column (time) has been resampled based on the selected method.

Examples

set.seed(123)
X <- cbind(runif(100), runif(100), sort(runif(100)))

# Pure permutation
sims_pure <- sim.procedures(X, nperm = 10, method = "pure")
head(sims_pure[[1]])
# Block permutation
sims_block <- sim.procedures(X, nperm = 10, nblocks = 5, method = "block")

# Visualize the first result from block permutation
plot_stpp(sims_block[[1]], type = "3D")