Help for package sppt

Type:

Package

Title:

Spatial Point Pattern Test for Areal Data

Version:

1.0.0

Description:

Implements the area-based nonparametric spatial point pattern test of Andresen (2009) <doi:10.1016/j.apgeog.2008.12.004> for comparing two spatial distributions of events aggregated to areal units, with the applications and extensions reviewed in Andresen (2016) <doi:10.1177/2059799116630659>. Unit-level inference uses bounded Monte Carlo bootstrap p-values, with marginal confidence interval overlap retained as an option. Global similarity is summarised by the robust S index, accompanied by a test of whether two patterns are more different than two realisations of a single process, and a test of the spatial arrangement of the differences. A longitudinal extension compares each period of a time series against a base period and classifies unit-level trajectories. Publication-ready maps can be redrawn with user-chosen layouts.

License:

MIT + file LICENSE

Encoding:

UTF-8

LazyData:

true

LazyDataCompression:

Depends:

R (≥ 4.1.0)

Imports:

classInt, ggplot2 (≥ 3.5.0), ggspatial, graphics, grDevices, Matrix, rlang, scales, sf, spdep, stats, tools, utils

Suggests:

knitr, rmarkdown, testthat (≥ 3.0.0)

VignetteBuilder:

knitr

Config/testthat/edition:

RoxygenNote:

7.3.3

NeedsCompilation:

Packaged:

2026-10-16 07:21:59 UTC; andre

Author:

Martin A. Andresen

[aut, cre, cph]

Maintainer:

Martin A. Andresen <andresen@sfu.ca>

Repository:

CRAN

Date/Publication:

2026-10-27 15:00:11 UTC

sppt: Spatial Point Pattern Test for Areal Data

Description

Compares spatial distributions of events aggregated to areal units.

Details

sppt() compares two patterns (a base and a test) unit by unit, using bounded Monte Carlo bootstrap p-values and, optionally, marginal confidence interval overlap.
sppt_longitudinal() compares every period of a time series against the first period and classifies unit-level trajectories.
sppt_sameness() tests a robust S index against sameness.
sppt_map() redraws the maps from either analysis with a chosen layout.

Author(s)

Maintainer: Martin A. Andresen andresen@sfu.ca (ORCID) [copyright holder]

References

Andresen, M.A. (2009). Testing for similarity in area-based spatial patterns: A nonparametric Monte Carlo approach. Applied Geography, 29(3), 333-345. doi:10.1016/j.apgeog.2008.12.004

Andresen, M.A. (2016). An area-based nonparametric spatial point pattern test: The test, its applications, and the future. Methodological Innovations, 9, Article 12. doi:10.1177/2059799116630659

Andresen, M.A. (2026). Testing longitudinal change in spatial distributions and patterns: An extension of the spatial point pattern test. Annals of the American Association of Geographers, in press.

Spatial point pattern test with bootstrap p-values

Description

Compares two spatial point patterns, a base and a test, across a set of areal units. Each unit is tested with a bounded Monte Carlo bootstrap p-value and, optionally, with marginal confidence interval overlap (Andresen 2009). Both are derived from a single shared set of B bootstrap replicates.

Global similarity is summarised by the robust S index, the proportion of eligible units (those not empty in both patterns) that are not flagged as different. It is accompanied by two further tests:

a sameness test (sppt_sameness()), which asks whether the two patterns are more different than two realisations of a single process;
an arrangement test, global Moran's I on the signed evidence surface \mathrm{sign}(d) \times -\log_{10}(p), which asks whether the differences, however many there are, are spatially organised. Local Moran clusters are outlined on the p-value map.

Usage

sppt(
  data,
  group_col = "group",
  count_col,
  B = 4999L,
  seed = NULL,
  alpha = 0.05,
  calculate_pvalues = TRUE,
  sameness_tau = 0.05,
  arrangement = TRUE,
  nsim_moran = 999L,
  calculate_ci = FALSE,
  conf_level = 0.683,
  fix_base = FALSE,
  use_percentages = TRUE,
  verbose = TRUE,
  create_maps = TRUE,
  map_crs = NULL,
  export_maps = FALSE,
  export_dir = NULL,
  map_dpi = 300L,
  export_results = FALSE,
  export_format = "gpkg",
  export_results_dir = NULL
)

Arguments

data

An sf object with one row per spatial unit.

group_col

Character. Name of the column that uniquely identifies each spatial unit (for example "DAUID" or "CTUID").

count_col

Character vector of length 2 giving the base and test count columns, in that order (for example c("BNER2003", "BNER2025")).

B

Integer. Number of bootstrap replicates. Default 4999; an odd value avoids ties.

seed

Integer or NULL. Random seed for reproducibility.

alpha

Numeric in (0, 1). Significance level for flagging a unit. Default 0.05.

calculate_pvalues

Logical. Compute bootstrap p-values. Default TRUE.

sameness_tau

Numeric in (0, 1). Significance level for the sameness test on the robust S index. Separate from alpha, which governs the flagging of individual units. Default 0.05.

arrangement

Logical. Run the arrangement test (global and local Moran's I on the signed evidence surface). Requires polygon geometry. Default TRUE.

nsim_moran

Integer. Permutations for the global Moran's I test. Default 999.

calculate_ci

Logical. Also compute marginal bootstrap confidence intervals and the CI-based S index (Andresen 2009). Default FALSE.

conf_level

Numeric in (0, 1). Confidence level for the interval overlap criterion. Default 0.683 (about one standard deviation per tail), which matches alpha = 0.05 p-value performance. Only used when calculate_ci = TRUE.

fix_base

Logical. If TRUE, the base pattern is not bootstrapped and its observed proportions are treated as fixed. This gives an asymmetric test, intended only for comparison against a known reference distribution. Default FALSE.

use_percentages

Logical. Express proportions as percentages (0-100) rather than as proportions (0-1). Default TRUE.

verbose

Logical. Print the test summary to the console. Default TRUE.

create_maps

Logical. Build the maps and print them. The maps are always stored on the result (see Value) when this is TRUE. Default TRUE.

map_crs

CRS for the maps, passed to sf::st_transform(). NULL uses the CRS of data.

export_maps

Logical. Save the maps as PNG files in export_dir. Default FALSE.

export_dir

Character. Directory for exported maps and results. Required when export_maps or export_results is TRUE; created if it does not exist.

map_dpi

Integer. Resolution of exported maps. Default 300.

export_results

Logical. Write the result to export_results_dir. Default FALSE.

export_format

Character. Format for exported results: "gpkg" (default), "shp", "csv", "txt" or "rds". Only "rds" keeps the summary and the maps.

export_results_dir

Character or NULL. Directory for exported results. Defaults to export_dir.

Details

P-values. For each unit the observed difference in proportions (test minus base) is compared with the centred bootstrap distribution of that difference. The p-value is the bounded Monte Carlo form

p = \frac{1 + r}{B + 1},

where r counts centred bootstrap differences at least as large in absolute value as the observed difference. The smallest attainable p-value is therefore 1/(B + 1). Units with no events in either pattern are excluded from inference and receive NA.

No multiple testing correction is applied. Simulation evidence shows that corrections remove power at fine spatial scales without a meaningful benefit in false discovery rate; the sameness and arrangement tests address the global question instead.

Power. lambda_min, the mean number of events per unit in the rarer pattern, governs power. Below 1.5 a low-power warning is issued and below 0.5 a critical warning; at those resolutions, interpret groups of adjacent flagged units rather than isolated ones, or aggregate to coarser units.

Value

An sf object of class sppt_result containing the columns of data plus:

base_prop, test_prop: Observed share of base and test events in each unit.
difference: test_prop minus base_prop.
lambda_min: Mean events per unit for the rarer pattern (the same for every row).
low_power_warning: TRUE when lambda_min < 1.5.
pvalue: Bootstrap p-value; NA for units empty in both patterns.
sig_pvalue: 1 if pvalue <= alpha, otherwise 0.
sig_pvalue_dir: 1 = test > base, -1 = base > test, 0 = not significant.
signed_evidence: \mathrm{sign}(d) \times -\log_{10}(p).
moran_hh, moran_ll, moran_cluster: Local Moran clusters of the signed evidence surface (arrangement test only).
ci_lower, ci_upper, ci_lower_test, ci_upper_test: Marginal confidence bounds for base and test (CI only).
intervals_overlap: 1 = intervals overlap (similar), 0 = they do not (CI only).
sig_ci_dir: Direction of CI-based differences (CI only).

The attribute "sppt_summary" holds the robust S index, the sameness and arrangement tests and the settings; summary() prints it. The attribute "sppt_maps" holds the maps as a named list of ggplot objects (difference, pvalue and, with CI, ci); use sppt_map() to redraw them with a different layout.

References

Andresen, M.A. (2009). Testing for similarity in area-based spatial patterns: A nonparametric Monte Carlo approach. Applied Geography, 29(3), 333-345. doi:10.1016/j.apgeog.2008.12.004

Andresen, M.A. (2016). An area-based nonparametric spatial point pattern test: The test, its applications, and the future. Methodological Innovations, 9, Article 12. doi:10.1177/2059799116630659

Examples

res <- sppt(vancouver_ct, group_col = "id",
            count_col = c("BNER2003", "BNER2025"),
            B = 199, seed = 42, create_maps = FALSE)
summary(res)
table(res$sig_pvalue_dir)

Longitudinal spatial point pattern test

Description

Extends the spatial point pattern test to a time series of three or more periods (Andresen 2026). Every period is compared with the first period (the base), and each unit's trajectory is classified.

Two parallel sets of results are produced:

Envelope results. Bootstrap confidence intervals are computed for every unit and period. A unit is stable when its intervals overlap in all periods; the largest set of mutually overlapping periods gives a partial-stability measure. Trends are fitted only for units that are not stable.
P-value results (calculate_pvalues = TRUE). Each period is tested against the base with the bounded Monte Carlo bootstrap p-value p = (1 + r)/(B + 1) used by sppt(). A unit is stable when no period is significant at alpha; trends are fitted only for units with at least one significant period.

Trends are ordinary least squares fits of each unit's share (or count) on time, linear or quadratic, tested at trend_alpha. Units are also grouped into Jenks natural-breaks classes of their base-period count, and trend outcomes are tabulated by class. Units that had events and then recorded none in each of the last min_zero_streak periods are flagged as sustained zeros.

Usage

sppt_longitudinal(
  data,
  group_col = "group",
  count_col,
  B = 4999L,
  seed = NULL,
  alpha = 0.05,
  calculate_pvalues = TRUE,
  conf_level = 0.683,
  fix_base = FALSE,
  use_percentages = TRUE,
  prefix = NULL,
  min_zero_streak = 3L,
  trend_type = "linear",
  trend_alpha = 0.05,
  jenks = "optimal",
  jenks_gvf_threshold = 0.05,
  verbose = TRUE,
  create_maps = TRUE,
  create_time_plot = TRUE,
  export_dir = NULL,
  export_maps = FALSE,
  map_format = "png",
  map_dpi = 300L,
  map_width_in = 7,
  map_height_in = 5.5,
  include_title = TRUE,
  figure_prefix = NULL,
  export_time_plot = FALSE,
  export_tables = FALSE,
  export_results = FALSE,
  export_format = "gpkg"
)

Arguments

data

An sf object with one row per spatial unit.

group_col

Character. Name of the column that uniquely identifies each spatial unit (for example "DAUID" or "CTUID").

count_col

Character vector of three or more count columns, one per period, in time order (for example paste0("BNER", 2003:2025)). The first is the base period. Column names should end in the period number, which is used to label the time axis.

B

Integer. Number of bootstrap replicates. Default 4999; an odd value avoids ties.

seed

Integer or NULL. Random seed for reproducibility.

alpha

Numeric in (0, 1). Significance level for flagging a unit. Default 0.05.

calculate_pvalues

Logical. Compute bootstrap p-values. Default TRUE.

conf_level

Numeric in (0, 1). Confidence level of the bootstrap envelopes. Default 0.683 (about one standard deviation per tail), which matches alpha = 0.05 p-value performance.

fix_base

Logical. If TRUE, the base period is not bootstrapped and its observed shares are treated as fixed (the approach of earlier longitudinal applications). Not recommended for standard use. Default FALSE.

use_percentages

Logical. Express proportions as percentages (0-100) rather than as proportions (0-1). Default TRUE.

prefix

Character or NULL. Prefix for the output column names, so that results for several crime types can sit in one object. NULL (default) derives it from count_col, for example "BNER_".

min_zero_streak

Integer. Number of final periods with no events that, after earlier events, mark a unit as a sustained zero. Default 3.

trend_type

Character. "linear" (default) or "quadratic". A quadratic trend needs four or more periods.

trend_alpha

Numeric in (0, 1). Significance level for the trend regressions; independent of alpha. Default 0.05.

jenks

"optimal" (default) to choose the number of Jenks classes by the goodness of variance fit (GVF) stopping rule, or an integer of at least 2 to fix it.

jenks_gvf_threshold

Numeric in (0, 1). With jenks = "optimal", classes are added until the GVF improves by less than this. Default 0.05.

verbose

Logical. Print the test summary to the console. Default TRUE.

create_maps

Logical. Build the maps and print them. Default TRUE.

create_time_plot

Logical. Plot total events per period with the summed envelope bounds. Default TRUE.

export_dir

Character. Directory for every export option. Required when any ⁠export_*⁠ argument is TRUE; created if it does not exist.

export_maps

Logical. Save the maps in export_dir. Default FALSE.

map_format

Character. File type for exported maps: "png" (default), "tif", "eps", "svg" or "pdf".

map_dpi

Integer. Resolution of exported maps and the time plot. Default 300; use 600 for journal TIFF submission.

map_width_in, map_height_in

Numeric. Size of exported maps in inches. Defaults 7 and 5.5.

include_title

Logical. Keep map titles in exported maps. FALSE removes them, as journals usually require. Default TRUE.

figure_prefix

Character or NULL. If given, exported maps are named ⁠<figure_prefix>1⁠, ⁠<figure_prefix>2⁠, and so on in map order (for example "SmithFigure") instead of descriptive names.

export_time_plot

Logical. Save the time plot as a PNG in export_dir. Default FALSE.

export_tables

Logical. Save the summary indices and the Jenks trend tables as a CSV in export_dir. Default FALSE.

export_results

Logical. Write the result to export_dir. Default FALSE.

export_format

Character. Format for exported results: "gpkg" (default), "shp", "csv", "txt" or "rds". Only "rds" keeps the summary and the maps.

Value

An sf object containing the columns of data plus, for each period, ⁠<count_col>_L⁠ and ⁠<count_col>_U⁠ (envelope bounds) and, with p-values, ⁠pvalue_<count_col>⁠ (NA for the base period). Unit-level results are named with prefix:

intervals_overlap: 1 if all envelopes overlap (stable).
max_overlap_proportion: Largest share of periods whose envelopes mutually overlap.
sustained_zero: 1 if the unit declined to zero and stayed there.
Trend: Trend class. Linear: 2 / -2 significant increase / decrease, 1 / -1 non-significant, 0 stable. Quadratic: 1 / -1 significant linear only, 2 / 3 significant increase with negative / positive curvature, -2 / -3 significant decrease with positive / negative curvature, 0 otherwise.
Trend_All, Trend_Linear_Sig: Linear direction, and significant linear direction (1, -1, 0).
Linear_R2, Linear_b0, Linear_bt, Linear_bt_sig: Linear fit; Linear_bt_sig is the slope for significant units only.
Quad_R2, Quad_b0, Quad_bt, Quad_bt2, Quad_Sig: Quadratic fit and the sign of a significant quadratic term (quadratic only).
⁠*_p⁠ versions: The same quantities gated on p-value significance: intervals_overlap_p, max_overlap_proportion_p (share of non-significant periods), prop_sig_p, sustained_zero_p, Trend_p, Trend_p_All, Trend_p_Linear_Sig, Linear_p_R2, Linear_p_b0, Linear_p_bt, Linear_p_bt_sig, Quad_p_R2, Quad_p_b0, Quad_p_bt, Quad_p_bt2, Quad_p_Sig.
jenks_level, jenks_min, jenks_max, baseline_interval: Jenks class of the base-period count, its bounds and label.

The S indices, Jenks breaks, trend tables and p-value matrix are stored as attributes (for example attr(x, "s_index_overlap"), attr(x, "interval_table_linear"), attr(x, "pvalue_matrix")). With create_maps = TRUE the maps are stored as a named list of ggplot objects in attr(x, "sppt_maps"); use sppt_map() to redraw them with a different layout.

References

Examples


res <- sppt_longitudinal(vancouver_ct, group_col = "id",
                         count_col = paste0("BNER", 2003:2025),
                         B = 199, seed = 42, verbose = FALSE,
                         create_maps = FALSE, create_time_plot = FALSE)
table(res$BNER_Trend_p)
attr(res, "si_overlap_p")

Redraw SPPT maps with a different layout

Description

sppt() and sppt_longitudinal() draw their maps with a fixed layout that suits Vancouver: north arrow top right, scale bar along the bottom, and legends placed where the city's shape leaves room. sppt_map() rebuilds any of those maps from a result object with the legend, north arrow and scale bar placed wherever suits your study area.

Every map is an ordinary ggplot object, so anything not covered by the arguments here can be changed by adding layers or a theme, for example sppt_map(res, "pvalue") + ggplot2::theme(legend.text = ggplot2::element_text(size = 12)).

Usage

sppt_map(
  result,
  which = NULL,
  legend_position = NULL,
  north_arrow = NULL,
  scale_bar = NULL,
  crs = NULL,
  title = TRUE
)

Arguments

result

A result returned by sppt() or sppt_longitudinal().

which

Character vector of map names. NULL (default) returns every map available for the result. For sppt() results: "difference", "pvalue" and "ci". For sppt_longitudinal() results: "overlap", "s_index_average", "sustained_zero", "trend_all", "trend_linear_sig", "linear_r2", "linear_slope", "quad_r2", "quad_trend", and the p-value versions "p_overlap", "p_sig_proportion", "p_sustained_zero", "p_trend_all", "p_trend_linear_sig", "p_linear_r2", "p_linear_slope", "p_quad_r2", "p_trend". The quadratic maps exist only when trend_type = "quadratic".

legend_position

Where to put the legend: "right", "left", "top", "bottom", "none", or a numeric vector c(x, y) giving the top-left corner of the legend in panel coordinates (0 to 1; negative x pushes it into the left margin). NULL (default) keeps each map's built-in placement.

north_arrow

Location of the north arrow: "tr", "tl", "br", "bl", or "none" to omit it. NULL (default) uses "tr".

scale_bar

Location of the scale bar, with the same options as north_arrow. NULL (default) uses "bl" for sppt() maps and "br" for sppt_longitudinal() maps.

crs

CRS to project the maps to, passed to sf::st_transform(). NULL keeps the CRS of the result.

title

Logical. Keep the map titles. FALSE removes them, as journals usually require. Default TRUE.

Value

A ggplot object when which names a single map, otherwise a named list of ggplot objects.

Examples

res <- sppt(vancouver_ct, group_col = "id",
            count_col = c("BNER2003", "BNER2025"),
            B = 199, seed = 42, verbose = FALSE, create_maps = FALSE)
p <- sppt_map(res, "pvalue", legend_position = "right",
              north_arrow = "tl", scale_bar = "br")
p

Test a robust S index against sameness

Description

Tests whether two spatial patterns are more different than two realisations of a single process would be. The null distribution of the S index follows from the calibration of the test itself: with units flagged at alpha, a proportion alpha of eligible units flag by chance under identity, so the number flagged is \mathrm{Binomial}(n, \alpha) and S is centred at 1 - \alpha with standard error \sqrt{\alpha(1 - \alpha)/n}, where n is the number of eligible units.

Verified against 10,120 simulated combinations (8 crime types, 23 years, 55 spatial scales; 5,000 iterations each). The mean error of the closed-form threshold is +0.0099 at alpha = 0.05 across the full grid, and +0.0006 where at least 95 percent of units are eligible. Where many units are empty in both patterns the test is mildly conservative, because discreteness in sparse units lowers the realised flagging rate below nominal.

sppt() calls this function automatically; it is exported so that the test can be applied to S indices computed elsewhere.

Usage

sppt_sameness(
  s_index,
  n_sig,
  n_eligible,
  alpha = 0.05,
  tau = 0.05,
  exact = NULL
)

Arguments

s_index

Numeric. Observed robust S index.

n_sig

Integer. Number of units flagged as dissimilar. Used for the exact binomial tail; may be NA, in which case the normal approximation is used.

n_eligible

Integer. Number of units in the S index denominator, that is, the total less units empty in both patterns.

alpha

Numeric in (0, 1). The unit-level significance level used in the spatial point pattern test. Default 0.05.

tau

Numeric in (0, 1). Significance level for the sameness test itself. This is a separate choice from alpha: alpha governs when an individual unit is flagged, tau governs how much of the lower tail of the null distribution of S is treated as ordinary. Default 0.05.

exact

Logical or NULL. Use the exact binomial tail rather than the normal approximation. NULL (default) selects the exact test when n_eligible < 200, where the two diverge materially and the normal approximation is anti-conservative.

Value

A list with elements s_index, n_eligible, alpha, tau, s_expected, se, threshold, z, p_value, method, and same (TRUE when the two patterns are no more different than two realisations of one process).

Examples

sppt_sameness(s_index = 0.90, n_sig = 20, n_eligible = 200, alpha = 0.05)

Summarise an SPPT result

Description

Prints the robust S index, the sameness and arrangement tests and the settings of a sppt() result.

Usage

## S3 method for class 'sppt_result'
summary(object, ...)

## S3 method for class 'sppt_result'
print(x, ...)

Arguments

object, x

An object of class sppt_result, as returned by sppt().

...

Further arguments passed to the sf method.

Value

summary() invisibly returns the summary list stored in the "sppt_summary" attribute. print() invisibly returns x.

Examples

res <- sppt(vancouver_ct, group_col = "id",
            count_col = c("BNER2003", "BNER2025"),
            B = 199, seed = 1, verbose = FALSE, create_maps = FALSE)
summary(res)

Residential and commercial burglary in Vancouver, 2003-2025

Description

Annual counts of residential burglary (break and enter, residential) and commercial burglary (break and enter, commercial) in the City of Vancouver, British Columbia, from 2003 to 2025, aggregated to 2021 Census spatial units.

vancouver_ct: 127 census tracts.
vancouver_da: 1,016 dissemination areas.

Census tracts are coarse enough that most period comparisons have good power, which makes vancouver_ct the better choice for examples and first analyses. Dissemination areas show the fine-scale behaviour of the test, including the low-power warnings.

Usage

vancouver_ct

vancouver_da

Format

An sf object with polygon geometry in NAD83(CSRS) / UTM zone 10N and the following columns:

id: Integer. Sequential unit identifier, for use as group_col. It is not a Statistics Canada code.
BNER2003, ..., BNER2025: Integer. Residential burglaries in each year.
BNEC2003, ..., BNEC2025: Integer. Commercial burglaries in each year.
geometry: sfc_MULTIPOLYGON boundaries of the unit.

vancouver_ct has 127 rows and vancouver_da has 1,016 rows.

An object of class sf (inherits from data.frame) with 1016 rows and 48 columns.

Source

Crime counts: Vancouver Police Department open crime data, https://geodash.vpd.ca/opendata/. Boundaries: Statistics Canada, 2021 Census boundary files, https://www12.statcan.gc.ca/census-recensement/2021/geo/sip-pis/boundary-limites/index2021-eng.cfm. Events were aggregated to the 2021 units by the package author.

Examples

vancouver_ct
colSums(sf::st_drop_geometry(vancouver_ct)[, c("BNER2003", "BNER2025")])
plot(vancouver_ct["BNER2025"])

Package {sppt}

sppt: Spatial Point Pattern Test for Areal Data

Description

Details

Author(s)

References

Spatial point pattern test with bootstrap p-values

Description

Usage

Arguments

Details

Value

References

See Also

Examples

Longitudinal spatial point pattern test

Description

Usage

Arguments

Value

References

See Also

Examples

Redraw SPPT maps with a different layout

Description

Usage

Arguments

Value

See Also

Examples

Test a robust S index against sameness

Description

Usage

Arguments

Value

See Also

Examples

Summarise an SPPT result

Description

Usage

Arguments

Value

Examples

Residential and commercial burglary in Vancouver, 2003-2025

Description

Usage

Format

Source

Examples