Package {spacc}

Title:	Fast Spatial Species Accumulation Curves
Version:	0.11.0
Copyright:	See inst/COPYRIGHTS for the copyright holders of the bundled nanoflann library.
Description:	High-performance spatial species accumulation curves using fixed-focus spatially constrained rarefaction, nearest-centroid expansion, and nearest-neighbor walks with a C++ backend for speed. Supports Hill numbers (q=0,1,2), beta diversity partitioning (turnover/nestedness) at site and pooled resolution, coverage-based rarefaction and extrapolation of alpha and beta diversity, phylogenetic diversity (PD; Faith's PD, EvoHeritage, mean pairwise distance, mean nearest taxon distance), functional diversity accumulation including aggregated functional diversity and intraspecific trait variation, phylogenetic and functional redundancy, interaction network accumulation, diversity-area relationships (DAR), endemism-area curves, sampling-effort correction and fragmentation analysis, species-area relationship (SAR) models based on extreme value theory (EVT), and area-based richness extrapolation via the total-species curve. Total richness is also estimated from frequency counts, including the accumulation rate curve estimator and an estimator for sites sampled without replacement, and shift-and-rotate null models test the position of a sampling design. Multiple starting points (seeds) provide uncertainty quantification. Methods are described in Chao et al. (2014) <doi:10.1890/13-0133.1>, Baselga (2010) <doi:10.1111/j.1466-8238.2009.00490.x>, Chao and Jost (2012) <doi:10.1890/11-1952.1>, Faith (1992) <doi:10.1016/0006-3207(92)91201-3>, Ma (2018) <doi:10.1002/ece3.4425>, Borda-de-Agua et al. (2025) <doi:10.1038/s41467-025-59239-7>, Hanski et al. (2013) <doi:10.1073/pnas.1311190110>, Jost (2007) <doi:10.1890/06-1736.1>, Ugland et al. (2003) <doi:10.1046/j.1365-2656.2003.00748.x>, Chao et al. (2023) <doi:10.1002/ecm.1588>, Chiu (2023) <doi:10.1111/2041-210X.14146>, Chiu et al. (2023) <doi:10.1098/rstb.2022.0183>, Rosindell et al. (2024) <doi:10.1093/sysbio/syad072>, Wojcik et al. (2024) <doi:10.1111/2041-210X.14470>, Ridder et al. (2024) <doi:10.1111/2041-210X.14443>, Jones et al. (2025) <doi:10.1111/ddi.70080>, and Shestopaloff (2025) <doi:10.1007/s10651-025-00644-y>.
License:	MIT + file LICENSE
Language:	en-US
Encoding:	UTF-8
Depends:	R (≥ 3.5)
LinkingTo:	Rcpp, RcppParallel
Imports:	Rcpp, RcppParallel, stats, parallel
Suggests:	testthat (≥ 3.0.0), knitr, rmarkdown, svglite, ggplot2, cli, sf, areaOfEffect, ape, iNEXT, iNEXT.3D, iNEXT.beta3D, betapart, cluster, vegan
Config/testthat/edition:	3
VignetteBuilder:	knitr
SystemRequirements:	GNU make
URL:	https://gillescolling.com/spacc/, https://github.com/gcol33/spacc
BugReports:	https://github.com/gcol33/spacc/issues
Config/roxygen2/version:	8.1.0
NeedsCompilation:	yes
Packaged:	2026-10-10 20:27:46 UTC; Gilles Colling
Author:	Gilles Colling [aut, cre, cph], Marius Muja [ctb, cph] (nanoflann library (inst/include/nanoflann.hpp, BSD license)), David G. Lowe [ctb, cph] (nanoflann library (inst/include/nanoflann.hpp, BSD license)), Jose Luis Blanco [ctb, cph] (nanoflann library (inst/include/nanoflann.hpp, BSD license))
Maintainer:	Gilles Colling <gilles.colling051@gmail.com>
Repository:	CRAN
Date/Publication:	2026-10-11 13:40:12 UTC

spacc: Fast Spatial Species Accumulation Curves

Description

High-performance spatial species accumulation curves using fixed-focus spatially constrained rarefaction, nearest-centroid expansion, and nearest-neighbor walks with a C++ backend for speed. Supports Hill numbers (q=0,1,2), beta diversity partitioning (turnover/nestedness) at site and pooled resolution, coverage-based rarefaction and extrapolation of alpha and beta diversity, phylogenetic diversity (PD; Faith's PD, EvoHeritage, mean pairwise distance, mean nearest taxon distance), functional diversity accumulation including aggregated functional diversity and intraspecific trait variation, phylogenetic and functional redundancy, interaction network accumulation, diversity-area relationships (DAR), endemism-area curves, sampling-effort correction and fragmentation analysis, species-area relationship (SAR) models based on extreme value theory (EVT), and area-based richness extrapolation via the total-species curve. Total richness is also estimated from frequency counts, including the accumulation rate curve estimator and an estimator for sites sampled without replacement, and shift-and-rotate null models test the position of a sampling design. Multiple starting points (seeds) provide uncertainty quantification. Methods are described in Chao et al. (2014) doi:10.1890/13-0133.1, Baselga (2010) doi:10.1111/j.1466-8238.2009.00490.x, Chao and Jost (2012) doi:10.1890/11-1952.1, Faith (1992) doi:10.1016/0006-3207(92)91201-3, Ma (2018) doi:10.1002/ece3.4425, Borda-de-Agua et al. (2025) doi:10.1038/s41467-025-59239-7, Hanski et al. (2013) doi:10.1073/pnas.1311190110, Jost (2007) doi:10.1890/06-1736.1, Ugland et al. (2003) doi:10.1046/j.1365-2656.2003.00748.x, Chao et al. (2023) doi:10.1002/ecm.1588, Chiu (2023) doi:10.1111/2041-210X.14146, Chiu et al. (2023) doi:10.1098/rstb.2022.0183, Rosindell et al. (2024) doi:10.1093/sysbio/syad072, Wojcik et al. (2024) doi:10.1111/2041-210X.14470, Ridder et al. (2024) doi:10.1111/2041-210X.14443, Jones et al. (2025) doi:10.1111/ddi.70080, and Shestopaloff (2025) doi:10.1007/s10651-025-00644-y.

Author(s)

Maintainer: Gilles Colling gilles.colling051@gmail.com (ORCID) [copyright holder]

Authors:

• Gilles Colling gilles.colling051@gmail.com (ORCID) [copyright holder]

Other contributors:

• Marius Muja (nanoflann library (inst/include/nanoflann.hpp, BSD license)) [contributor, copyright holder]

• David G. Lowe (nanoflann library (inst/include/nanoflann.hpp, BSD license)) [contributor, copyright holder]

• Jose Luis Blanco (nanoflann library (inst/include/nanoflann.hpp, BSD license)) [contributor, copyright holder]

See Also

Useful links:

• https://gillescolling.com/spacc/

• https://github.com/gcol33/spacc

• Report bugs at https://github.com/gcol33/spacc/issues

ACE Richness Estimator

Description

Abundance-based Coverage Estimator (ACE) of total species richness (Chao & Lee 1992). Separates species into rare and abundant groups based on an abundance threshold.

Usage

ace(x, threshold = 10L)

Arguments

Details

ACE partitions species into rare (abundance \le threshold) and abundant (abundance > threshold) groups. The estimate is:

S_{ACE} = S_{abund} + \frac{S_{rare}}{C_{ACE}} + \frac{f_1}{C_{ACE}} \gamma^2

where C_{ACE} is the sample coverage of rare species and \gamma^2 is the estimated coefficient of variation.

Value

An object of class spacc_estimate.

References

Chao, A. & Lee, S.M. (1992). Estimating the number of classes via sample coverage. Journal of the American Statistical Association, 87, 210-217.

See Also

chao1(), chao2()

Examples

species <- matrix(rpois(50 * 30, 2), nrow = 50)
ace(species)

Alpha Diversity (Per-Site)

Description

Compute Hill numbers for each site individually.

Usage

alphaDiversity(x, q = c(0, 1, 2), coords = NULL)

Arguments

Details

Alpha diversity represents local (within-site) diversity. For Hill numbers:

• q=0: Species richness

• q=1: Exponential of Shannon entropy

• q=2: Inverse Simpson concentration

Value

If coords is NULL, a data.frame with columns for each q value. If coords is provided, a spacc_alpha object.

References

Jost, L. (2007). Partitioning diversity into independent alpha and beta components. Ecology, 88, 2427-2439.

See Also

gammaDiversity() for regional diversity, diversityPartition() for full alpha-beta-gamma decomposition

Examples

species <- matrix(rpois(50 * 30, 2), nrow = 50)
alpha <- alphaDiversity(species, q = c(0, 1, 2))
head(alpha)

# Mean alpha diversity
colMeans(alpha)

Analytical Accumulation Methods

Description

These methods compute expected species accumulation without simulation. They are faster but don't provide spatial information.

Accumulation Rate Curve (ARC) Richness Estimator

Description

Estimate total species richness from the rate at which new species appear along a random reordering of the observations (Shestopaloff 2025). The rate curve is extrapolated with three decay curves whose tails give the number of species not yet seen, combined with bounded weights.

Usage

arc(
  x,
  data = c("abundance", "incidence"),
  K = 20L,
  adjust = FALSE,
  n_boot = 2000L,
  parallel = TRUE,
  n_cores = NULL
)

Arguments

Details

Let y_k be the number of species observed k times and N the number of observations of the species with k \le K. Drawing the observations in random order, the expected number of these species seen after n draws is

S(n) = \sum_{k=1}^{K} y_k \left[1 - (1 - n/N)^k\right],

and the accumulation rate is s_n = S(n) - S(n-1), evaluated at up to 1000 points. Three curves are fitted to the tail of s_n by least squares with weights proportional to the squared point index: the log-normal 1 - \Phi((\log n - \mu)/\sigma), the log-logistic 1/(1 + \alpha n^\gamma) and the Weibull \exp(-\lambda n^\gamma), the last on the log scale. Each curve's integral from N to infinity is its number of undetected species; the log-logistic integral is finite only for \gamma > 1.

The estimate combines the three with weights fitted by the squared misfit to s_n plus a penalty on the product of the log-logistic and Weibull weights (Eq. 8), each bounded by caps built from the three estimates and the observed richness (Eq. 7). The undetected species are averaged over fits to the last 100, 95, ..., 70 percent of the curve.

The penalty is there to make the fit prefer either the log-logistic or the Weibull, and gives the loss one local minimum for each. The weights are the end point of the projected gradient flow of the loss started at zero weights, which arc() solves in closed form face by face of the box of admissible weights. This is the minimum the author's reference code reaches (L-BFGS within an augmented Lagrangian, started at zero). In Poisson-Gamma samples of intermediate sparsity it is the log-logistic one, with median bias between -1 and +4 percent; the global minimum of the loss is the Weibull one there, with median bias between -9 and -17 percent (vignette("richness-estimation")).

With adjust = TRUE the counts y_1, \dots, y_{k_1} (the first run of positive counts) are mixed half and half with a log-linear fit \log y_k = b_0 + b_1 k, and fully replaced by it; the adjusted estimate is 2 S_{0.5} - S_1. The adjustment assumes the ratios (k+1) y_{k+1} / y_k are linear in k.

The interval comes from a parametric bootstrap of the counts. Their covariance is that of a mixture of 30 Poisson distributions whose weights are fitted to (y_1, \dots, y_K, y_{K+}) by penalised least squares (Eq. 12-15). The penalty is chosen by resampling species: 1500 resamples are ordered by their distance to the left-out species and split into deciles, and the first decile whose out-of-sample error is smallest at a non-zero penalty sets it. n_boot count vectors are drawn from a normal distribution with that covariance, rounded to non-negative integers and re-estimated; the interval is their 2.5 and 97.5 percent quantiles and the standard error their standard deviation. In Poisson-Gamma simulations the 95 percent interval covered the truth in 93 to 96 percent of samples for sparse data and about 89 percent for dense data, where the point estimate is biased upwards (vignette("richness-estimation")).

With data = "incidence" the observations are species-site occurrences and y_k becomes Q_k, the number of species found at k sites, the substitution that turns Chao1 into Chao2.

Value

An object of class spacc_estimate with estimator = "ARC". Its details hold undetected (undetected species from each curve) and weights (their weights), both averaged over the curve portions, K, data, penalty (the selected Poisson-mixture penalty) and boot (the bootstrap estimates).

References

Shestopaloff, K. (2025). An accumulation rate curve estimator for total species. Environmental and Ecological Statistics, 32, 311-345. doi:10.1007/s10651-025-00644-y

See Also

chao1(), chao2(), ace()

Examples

set.seed(1)
species <- matrix(rpois(40 * 300, rgamma(300, 0.4, 2)), nrow = 40, byrow = TRUE)
arc(species, n_boot = 0)

arc(species, n_boot = 200)

Convert spacc_metrics to sf

Description

Convert metrics to an sf object for spatial analysis and integration with the areaOfEffect package.

Usage

as_sf(x, crs = NULL)

Arguments

Value

An sf object with POINT geometries and metric columns.

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
metrics <- spaccMetrics(species, coords)
if (requireNamespace("sf", quietly = TRUE)) {
  metrics_sf <- as_sf(metrics)
}

Convert an Accumulation Curve from Another Package to a spacc Object

Description

Bring an existing species accumulation curve into the spacc object system so it can be summarised, plotted, extrapolated and compared with compare() alongside spatial curves.

Usage

as_spacc(x, ...)

## Default S3 method:
as_spacc(x, ...)

## S3 method for class 'specaccum'
as_spacc(x, ...)

## S3 method for class 'numeric'
as_spacc(x, ...)

## S3 method for class 'data.frame'
as_spacc(x, column = NULL, ...)

Arguments

Details

A specaccum object built with method = "random" stores one richness curve per permutation in its perm matrix. Those curves become the rows of the returned curves matrix and behave exactly like spacc seeds, so the across-curve confidence band and the compare() permutation, bootstrap, and AUC tests carry over unchanged. A specaccum object without a permutation matrix (for example method = "exact" or "collector"), a numeric vector and a data frame column each contribute a single mean curve, entry k being the value after k sites.

Rarefy::directionalSAC(community, gradient) with a site distance matrix as gradient orders the sites by distance from each plot in turn and averages the collector curves; this is spacc(x, coords, method = "knn", focal_points = coords), so the two curves agree up to the order of sites at tied distances.

Coordinates cannot be recovered from these objects, so the result carries no coords: it is a curve container for plotting and comparison, not a re-runnable spatial object.

Value

An object of class spacc with print(), summary(), plot(), as.data.frame(), and compare() support.

References

Thouverai, E., Pavoine, S., Tordoni, E., Chiarucci, A., Ricotta, C., Rocchini, D. & Bacaro, G. (2026). Rarefy: an R package for the calculation of taxonomic, functional and phylogenetic sample-based spatially explicit rarefaction curves. Ecological Informatics, 103803. doi:10.1016/j.ecoinf.2026.103803

See Also

compare() to test a converted curve against a spatial one.

Examples

if (requireNamespace("vegan", quietly = TRUE)) {
  data(BCI, package = "vegan")
  classic <- vegan::specaccum(BCI, method = "random")
  sp <- as_spacc(classic)
  summary(sp)
  plot(sp)
}

as_spacc(c(5, 8, 10, 11, 12))

Autoplot Methods for spacc Objects

Description

ggplot2::autoplot() methods for spacc objects. Each one calls the corresponding plot() method, which already returns a ggplot2 object.

Usage

## S3 method for class 'spacc'
autoplot(object, ...)

Arguments

Value

A ggplot2 object.

Examples

if (requireNamespace("ggplot2", quietly = TRUE)) {
  coords <- data.frame(x = runif(30), y = runif(30))
  species <- matrix(rbinom(30 * 15, 1, 0.3), nrow = 30)
  sac <- spacc(species, coords, n_seeds = 5, progress = FALSE)
  ggplot2::autoplot(sac)
}

Beta Distance-Decay

Description

Compute the distance-decay of community similarity by fitting models to pairwise beta dissimilarity as a function of geographic distance.

Usage

betaDecay(
  x,
  coords,
  index = c("sorensen", "jaccard"),
  model = c("all", "exponential", "power", "linear"),
  distance = c("euclidean", "haversine"),
  n_bins = 50L,
  progress = TRUE
)

Arguments

Details

Three decay models are available:

• Exponential: \beta = 1 - a \cdot e^{-b \cdot d}

• Power: \beta = a \cdot d^b

• Linear: \beta = a + b \cdot d

The half-life (exponential model) is the distance at which similarity decays to half its initial value: d_{1/2} = \ln(2) / b.

Value

An object of class spacc_beta_decay containing:

References

Nekola, J.C. & White, P.S. (1999). The distance decay of similarity in biogeography and ecology. Journal of Biogeography, 26, 867-878.

Morlon, H., Chuyong, G., Condit, R., et al. (2008). A general framework for the distance-decay of similarity in ecological communities. Ecology Letters, 11, 904-917.

See Also

spaccBeta() for spatial beta diversity accumulation, distanceDecay() for species similarity decay

Examples

coords <- data.frame(x = runif(30), y = runif(30))
species <- matrix(rbinom(30 * 20, 1, 0.3), nrow = 30)

bd <- betaDecay(species, coords)
print(bd)
plot(bd)

Bootstrap Richness Estimator

Description

Bootstrap estimator of total species richness (Smith & van Belle 1984). Uses species detection probabilities to estimate undetected species.

Usage

bootstrap_richness(x, n_boot = 200L)

Arguments

Details

The bootstrap estimator is:

S_{boot} = S_{obs} + \sum_{i=1}^{S_{obs}} (1 - p_i)^n

where p_i is the proportion of sites where species i occurs.

Value

An object of class spacc_estimate.

References

Smith, E.P. & van Belle, G. (1984). Nonparametric estimation of species richness. Biometrics, 40, 119-129.

See Also

jackknife(), chao2()

Examples

species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
bootstrap_richness(species, n_boot = 100)

Combine spacc Objects

Description

Combine spacc Objects

Usage

## S3 method for class 'spacc'
c(...)

Arguments

Value

A grouped spacc object with per-group curves.

Combine spacc_beta Objects

Description

Combine multiple spacc_beta objects by stacking their curve matrices. All objects must have the same number of sites and index.

Usage

## S3 method for class 'spacc_beta'
c(...)

Arguments

Value

A combined spacc_beta object with more seeds.

Combine spacc_coverage Objects

Description

Combine multiple spacc_coverage objects by stacking their curve matrices. All objects must have the same number of sites.

Usage

## S3 method for class 'spacc_coverage'
c(...)

Arguments

Value

A combined spacc_coverage object with more seeds.

Combine spacc_hill Objects

Description

Combine multiple spacc_hill objects by stacking their curve matrices. All objects must have the same number of sites.

Usage

## S3 method for class 'spacc_hill'
c(...)

Arguments

Value

A combined spacc_hill object with more seeds.

Chao1 Richness Estimator

Description

Estimate total species richness from abundance data using the Chao1 estimator (Chao 1984). Uses the number of singletons and doubletons to estimate undetected species.

Usage

chao1(x)

Arguments

Details

The Chao1 estimator is:

S_{Chao1} = S_{obs} + \frac{f_1^2}{2 f_2}

where f_1 is the number of singletons (species with total abundance 1) and f_2 is the number of doubletons (abundance 2).

When f_2 = 0, the bias-corrected form is used:

S_{Chao1} = S_{obs} + \frac{f_1 (f_1 - 1)}{2}

Value

An object of class spacc_estimate with components:

estimator

Name of the estimator ("chao1")

estimate

Estimated total richness

se

Standard error of the estimate

lower

Lower 95 percent confidence bound

upper

Upper 95 percent confidence bound

S_obs

Observed species richness

details

List with f1 (singletons) and f2 (doubletons)

References

Chao, A. (1984). Nonparametric estimation of the number of classes in a population. Scandinavian Journal of Statistics, 11, 265-270.

Chao, A. (1987). Estimating the population size for capture-recapture data with unequal catchability. Biometrics, 43, 783-791.

See Also

chao2() for incidence-based estimation, ace() for abundance-based coverage estimator

Examples

species <- matrix(rpois(50 * 30, 2), nrow = 50)
chao1(species)

Chao2 Richness Estimator

Description

Estimate total species richness from incidence (presence/absence) data using the Chao2 estimator (Chao 1987).

Usage

chao2(x, n_units = NULL)

Arguments

Details

The Chao2 estimator is the incidence-based analogue of Chao1:

S_{Chao2} = S_{obs} + \frac{t-1}{t} \frac{Q_1^2}{2 Q_2}

where t is the number of sites, Q_1 is the number of uniques (species found at exactly 1 site) and Q_2 is the number of duplicates (species found at exactly 2 sites).

Chao2 assumes the sites are drawn with replacement, or from an effectively infinite set of units. When they are a large fraction w = t / T of a finite set of T units it overestimates richness and does not reach the true richness as w \to 1. With n_units = T, chao2() returns the moment-based estimator for sampling without replacement (MoRE, Chiu 2023, Eq. 5):

S_{MoRE} = S_{obs} + \frac{t-1}{t} \frac{(1-w) Q_1^2}{Q_2 + w (Q_1 - Q_2)^+}

where (A)^+ is A for A > 0 and 1 otherwise; for Q_2 = 0 the undetected term is \frac{t-1}{t} \frac{1-w}{w} (Q_1 - 1). MoRE equals the observed richness at w = 1. Its standard error is the delta-method variance of Chiu (2023) and the interval is the log-normal one used for Chao2.

Value

An object of class spacc_estimate.

References

Chao, A. (1987). Estimating the population size for capture-recapture data with unequal catchability. Biometrics, 43, 783-791.

Chiu, C.-H. (2023). A species richness estimator for sample-based incidence data sampled without replacement. Methods in Ecology and Evolution, 14, 2482-2493. doi:10.1111/2041-210X.14146

See Also

chao1() for abundance-based estimation

Examples

species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
chao2(species)
# the 50 sites are 80% of a 62-cell grid
chao2(species, n_units = 62)

Coleman Expected Accumulation

Description

Compute the expected species accumulation curve using the Coleman method (Coleman et al. 1982). This is an analytical formula, no simulation needed.

Usage

coleman(x)

Arguments

Value

A data.frame with columns: sites, expected, sd

References

Coleman, B.D., Mares, M.A., Willig, M.R. & Hsieh, Y.H. (1982). Randomness, area, and species richness. Ecology, 63, 1121-1133.

Collector's Curve

Description

Compute the species accumulation curve in the order sites appear in the data (no randomization). Useful for understanding how data was collected.

Usage

collector(x)

Arguments

Value

A data.frame with columns: sites, species

Compare Two Accumulation Curves

Description

Test whether two species accumulation curves differ significantly.

Usage

compare(
  x,
  y,
  method = c("permutation", "bootstrap", "auc"),
  normalize = FALSE,
  n_perm = 999L,
  ...
)

Arguments

Value

An object of class spacc_comp containing:

References

Colwell, R.K., Mao, C.X. & Chang, J. (2004). Interpolating, extrapolating, and comparing incidence-based species accumulation curves. Ecology, 85, 2717-2727.

Gotelli, N.J. & Colwell, R.K. (2011). Estimating species richness. In: Biological Diversity: Frontiers in Measurement and Assessment, pp. 39-54. Oxford University Press.

Examples

coords <- data.frame(x = runif(50), y = runif(50))
sp_a <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
sp_b <- matrix(rbinom(50 * 30, 1, 0.5), nrow = 50)

sac_a <- spacc(sp_a, coords, n_seeds = 10)
sac_b <- spacc(sp_b, coords, n_seeds = 10)
comp <- compare(sac_a, sac_b)
print(comp)

Compare Multiple SAR Models

Description

Fit all (or a subset of) asymptotic species-area models and compare them using AIC, BIC, delta-AIC, and Akaike weights.

Usage

compareModels(
  object,
  models = c("michaelis-menten", "lomolino", "asymptotic", "weibull", "logistic", "evt"),
  ...
)

Arguments

Value

An object of class spacc_model_compare containing:

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
sac <- spacc(species, coords)
cm <- compareModels(sac)
print(cm)

Sample Completeness Profile

Description

Compute the ratio of observed to estimated diversity across diversity orders, measuring how complete a sample is at each level of the Hill number spectrum.

Usage

completenessProfile(x, q = seq(0, 2, by = 0.2), coords = NULL)

Arguments

Details

Sample completeness at order q is:

C_q = \frac{D_q^{obs}}{D_q^{est}}

where D_q^{obs} is the observed Hill number and D_q^{est} is the estimated asymptotic Hill number.

Completeness near 1 means the sample captures most of the true diversity at that order. Completeness typically increases with q because dominant species are detected early.

Asymptotic estimators used:

• q = 0: Chao1 estimator

• q = 1: Chao & Jost (2015) entropy estimator

• q = 2: Inverse Simpson estimator with bias correction

• Other q: Interpolated between adjacent integer estimates

When coords is provided, per-site completeness is computed by treating each site's abundance vector as an independent sample.

Value

An object of class spacc_completeness containing:

References

Chao, A. & Jost, L. (2012). Coverage-based rarefaction and extrapolation: standardizing samples by completeness rather than size. Ecology, 93, 2533-2547.

Chao, A. & Jost, L. (2015). Estimating diversity and entropy profiles via discovery rates of new species. Methods in Ecology and Evolution, 6, 873-882.

See Also

diversityProfile() for observed diversity profiles, chao1() for richness estimation

Examples

species <- matrix(rpois(50 * 30, 2), nrow = 50)
comp <- completenessProfile(species)
print(comp)


plot(comp)

Confidence interval for an extrapolated asymptote

Description

Returns the calibrated bootstrap interval stored on the fit by default. With method = "profile" the (over-confident) nls profile interval is computed instead. The interval is for the asymptote parameter; see the bias caveat in extrapolate().

Usage

## S3 method for class 'spacc_fit'
confint(object, parm, level = 0.95, method = c("stored", "profile"), ...)

Arguments

Value

A 1-row matrix with the lower and upper bound for a.

Diversity-Area Relationship (DAR)

Description

Diversity against cumulative area along a spatial accumulation, for Hill numbers (Ma 2018) or for any curve spacc computes: species richness, phylogenetic and functional metrics, redundancy.

Usage

dar(x, ...)

## Default S3 method:
dar(
  x,
  coords,
  q = c(0, 1, 2),
  n_seeds = 50L,
  method = "knn",
  area_method = c("convex_hull", "voronoi", "count"),
  distance = c("euclidean", "haversine"),
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  seed = NULL,
  focal_points = NULL,
  focal_domain = NULL,
  ...
)

## S3 method for class 'spacc'
dar(x, area_method = c("convex_hull", "voronoi", "count"), ...)

## S3 method for class 'spacc_hill'
dar(x, area_method = c("convex_hull", "voronoi", "count"), ...)

## S3 method for class 'spacc_phylo'
dar(x, area_method = c("convex_hull", "voronoi", "count"), ...)

## S3 method for class 'spacc_func'
dar(x, area_method = c("convex_hull", "voronoi", "count"), ...)

## S3 method for class 'spacc_redundancy'
dar(x, area_method = c("convex_hull", "voronoi", "count"), ...)

Arguments

Details

The DAR (Ma, 2018) generalizes the SAR by replacing species richness (q=0) with Hill numbers of any order. This reveals how different facets of diversity scale with area:

• q=0 (richness) recovers the classic SAR

• q=1 (Shannon) shows how common species diversity scales

• q=2 (Simpson) shows how dominant species diversity scales

The area is accumulated along the same site orders as the diversity curve, so any curve object gives its own diversity-area relationship: functional richness, dispersion and other facets can scale with area differently from richness (Dias et al. 2025). The power-law exponent z is the slope of log mean diversity against log mean area over the steps where both are positive.

Value

An object of class spacc_dar containing:

References

Ma, Z.S. (2018). DAR (diversity-area relationship): extending classic SAR for biodiversity and biogeography analyses. Ecology and Evolution, 8, 10023-10038.

Arrhenius, O. (1921). Species and area. Journal of Ecology, 9, 95-99.

Dias, R.A., Bastazini, V.A.G., Gianuca, A.T. & Matthews, T.J. (2025). The spatial scaling of multiple dimensions of functional diversity in habitat islands. Functional Ecology. doi:10.1111/1365-2435.70158

See Also

spaccHill(), spaccFunc(), spaccPhylo(), extrapolate()

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rpois(50 * 30, 2), nrow = 50)

d <- dar(species, coords, q = c(0, 1, 2))
plot(d)

traits <- matrix(rnorm(30 * 2), nrow = 30)
func <- spaccFunc(species, coords, traits, metric = c("fdis", "fric"))
dar(func)

Distance-Decay Analysis

Description

Analyze how species richness changes with distance from focal points.

Usage

distanceDecay(
  x,
  coords,
  n_seeds = 50L,
  breaks = NULL,
  distance = c("euclidean", "haversine"),
  progress = TRUE,
  seed = NULL
)

Arguments

Value

An object of class spacc_decay with distance-species relationships.

References

Nekola, J.C. & White, P.S. (1999). The distance decay of similarity in biogeography and ecology. Journal of Biogeography, 26, 867-878.

Soininen, J., McDonald, R. & Hillebrand, H. (2007). The distance decay of similarity in ecological communities. Ecography, 30, 3-12.

Compute Distance Matrix

Description

Pre-compute pairwise distances between sites for reuse across multiple spacc() calls. Supports sf objects with accurate geodesic distances for global-scale studies.

Usage

distances(x, method = NULL, fun = NULL, which = NULL)

Arguments

Details

For continental and global-scale studies, use sf objects with geographic CRS (e.g., EPSG:4326). The function will automatically use accurate geodesic distances via the S2 spherical geometry library.

For smaller study areas with projected coordinates (UTM, etc.), Euclidean distance is appropriate and faster.

Value

An object of class spacc_dist containing the distance matrix with coordinates stored as an attribute.

Examples

coords <- data.frame(x = runif(50), y = runif(50))
d <- distances(coords)

# Reuse for multiple analyses
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
sac <- spacc(species, d)

Alpha-Beta-Gamma Diversity Partitioning

Description

Decompose regional (gamma) diversity into local (alpha) and turnover (beta) components using multiplicative partitioning of Hill numbers.

Usage

diversityPartition(x, q = c(0, 1, 2), weights = "equal", coords = NULL)

Arguments

Details

This function implements Jost (2007) multiplicative partitioning:

\gamma = \alpha \times \beta

Where:

• Alpha: Mean effective number of species per site

• Beta: Effective number of distinct communities (1 = all identical, n_sites = all completely different)

• Gamma: Total effective number of species in the region

Beta diversity is interpreted as the effective number of communities:

• Beta = 1: All sites have identical composition

• Beta = n_sites: Sites share no species

Value

An object of class spacc_partition containing:

References

Jost, L. (2007). Partitioning diversity into independent alpha and beta components. Ecology, 88, 2427-2439.

Chao, A., Chiu, C.H. & Jost, L. (2014). Unifying species diversity, phylogenetic diversity, functional diversity, and related similarity and differentiation measures through Hill numbers. Annual Review of Ecology, Evolution, and Systematics, 45, 297-324.

See Also

alphaDiversity(), gammaDiversity(), spaccBeta() for spatial beta diversity accumulation

Examples

# Simulated community data
species <- matrix(rpois(50 * 30, 2), nrow = 50)

# Partition diversity
part <- diversityPartition(species, q = c(0, 1, 2))
print(part)

# Beta near 1 = homogeneous region
# Beta near n_sites = heterogeneous region

Diversity Profile

Description

Compute Hill numbers across a continuous range of diversity orders (q), producing a diversity profile for each site and the regional pool.

Usage

diversityProfile(
  x,
  q = seq(0, 3, by = 0.1),
  type = c("both", "per_site", "regional"),
  traits = NULL,
  tree = NULL,
  dist_method = c("euclidean", "gower"),
  normalize = TRUE,
  coords = NULL
)

Arguments

Details

A diversity profile plots effective number of species as a function of the order q. The key property is that Hill numbers are non-increasing in q: D_q \ge D_{q'} for q \le q'.

• q = 0: Species richness (insensitive to abundance)

• q = 1: Exponential of Shannon entropy (all species weighted equally by their proportional abundance)

• q = 2: Inverse Simpson concentration (emphasizes dominant species)

• q > 2: Increasingly dominated by common species

Supplying traits computes functional Hill numbers via a trait-similarity matrix (Leinster & Cobbold 2012); supplying tree computes phylogenetic Hill numbers weighted by branch length (Chao et al. 2010). Both return a spacc_profile whose profile_type records which was computed.

Value

An object of class spacc_profile containing:

References

Leinster, T. & Cobbold, C.A. (2012). Measuring diversity: the importance of species similarity. Ecology, 93, 477-489.

Chao, A., Chiu, C.H. & Jost, L. (2014). Unifying species diversity, phylogenetic diversity, functional diversity, and related similarity and differentiation measures through Hill numbers. Annual Review of Ecology, Evolution, and Systematics, 45, 297-324.

See Also

alphaDiversity() for per-site values at specific q, gammaDiversity() for regional diversity, evenness() for evenness profiles

Examples

species <- matrix(rpois(50 * 30, 2), nrow = 50)
prof <- diversityProfile(species)
print(prof)


plot(prof)

# Functional profile (pass a species-by-traits data.frame)
colnames(species) <- paste0("sp", 1:30)
traits <- data.frame(size = rnorm(30), row.names = paste0("sp", 1:30))
fp <- diversityProfile(species, traits = traits)

Evenness Profiles

Description

Compute species evenness across sites using Hill-based, Pielou, or Simpson evenness measures.

Usage

evenness(
  x,
  q = seq(0.1, 3, by = 0.1),
  type = c("hill", "pielou", "simpson"),
  coords = NULL
)

Arguments

Details

All evenness measures are bounded in [0, 1]:

• 0 = maximally uneven (one dominant species)

• 1 = perfectly even (all species equally abundant)

Hill evenness (Jost 2010):

E_q = D_q / D_0

This divides the effective number of species at order q by species richness.

Pielou's J (Pielou 1966):

J = \frac{\log(D_1)}{\log(S)} = \frac{H'}{\log(S)}

Simpson evenness:

E_{1/D} = \frac{1}{D_2 \cdot S}

Value

An object of class spacc_evenness containing:

References

Jost, L. (2010). The relation between evenness and diversity. Diversity, 2, 207-232.

Pielou, E.C. (1966). The measurement of diversity in different types of biological collections. Journal of Theoretical Biology, 13, 131-144.

See Also

diversityProfile() for Hill number profiles, alphaDiversity() for raw diversity values

Examples

species <- matrix(rpois(50 * 30, 2), nrow = 50)

# Hill evenness profile
ev <- evenness(species)
print(ev)

# Pielou's J
ev_j <- evenness(species, type = "pielou")
print(ev_j)

Extrapolate Total Species Richness

Description

Fit an asymptotic model to the spatial accumulation curve and estimate total species richness beyond the observed sampling effort.

Usage

extrapolate(
  object,
  model = c("michaelis-menten", "lomolino", "asymptotic", "weibull", "logistic", "evt"),
  interval = c("bootstrap", "profile", "none"),
  R = 200L,
  level = 0.95,
  compare = TRUE,
  n_units = NULL,
  warn_ratio = 2,
  ...
)

Arguments

Value

An object of class spacc_fit containing:

Bias caveat

The asymptote is a parametric extrapolation of the accumulation curve. On clustered or strongly under-sampled presence-absence data it tends to overestimate true richness, sometimes substantially, because the curve is far from saturation. The bootstrap interval quantifies the uncertainty of the fitted asymptote (curve-fit and across-seed variability) and is much wider than the nls profile interval, but it is not a calibrated interval for true total richness: it is centred on a possibly biased point estimate. For calibrated total-richness estimates prefer the nonparametric estimators chao2() / iChao2(), which are unbiased on the same data. The printout shows their values alongside the asymptote for comparison.

References

Lomolino, M.V. (2000). Ecology's most general, yet protean pattern: the species-area relationship. Journal of Biogeography, 27, 17-26.

Flather, C.H. (1996). Fitting species-accumulation functions and assessing regional land use impacts on avian diversity. Journal of Biogeography, 23, 155-168.

Borda-de-Agua, L., Whittaker, R.J., Cardoso, P., et al. (2025). Extreme value theory explains the topography and scaling of the species-area relationship. Nature Communications, 16, 5346.

See Also

extrapolateArea() for area-based extrapolation to a region larger than the one sampled; chao2(), iChao2() for nonparametric richness.

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
sac <- spacc(species, coords)
fit <- extrapolate(sac)
print(fit)

Extrapolate Richness to a Larger Area (Ugland T-S Curve)

Description

Estimate species richness for a spatial extent larger than the one sampled, using the total-species (T-S) curve of Ugland, Gray & Ellingsen (2003). Sites are partitioned into spatial subareas; the T-S curve is the expected total richness of random combinations of subareas plotted against their (convex-hull) area; an asymptotic model is fitted to that curve and extrapolated to a target area, with a site-bootstrap confidence interval.

Usage

extrapolateArea(
  x,
  coords,
  target_area = NULL,
  n_subareas = 10L,
  model = c("michaelis-menten", "asymptotic", "weibull", "lomolino", "logistic"),
  n_combos = 30L,
  R = 200L,
  level = 0.95,
  seed = NULL,
  progress = TRUE
)

Arguments

Details

Unlike extrapolate(), which extends the curve in sample count, this answers "how many species in a region of a given area", including areas beyond the sampled extent.

Value

An object of class spacc_area with components:

References

Ugland, K.I., Gray, J.S. & Ellingsen, K.E. (2003). The species-accumulation curve and estimation of species richness. Journal of Animal Ecology, 72, 888-897. doi:10.1046/j.1365-2656.2003.00748.x

See Also

extrapolate() for sample-count extrapolation; chao2() for nonparametric richness.

Examples

coords <- data.frame(x = runif(150), y = runif(150))
species <- matrix(rbinom(150 * 60, 1, 0.15), nrow = 150)
area_fit <- extrapolateArea(species, coords, n_subareas = 8, progress = FALSE)
print(area_fit)

Extrapolate Diversity Beyond Observed Coverage

Description

Hill numbers at target sample coverages along a coverage-tracking accumulation, interpolated within the observed coverage and extrapolated beyond it (Chao et al. 2014).

Usage

extrapolateCoverage(x, target_coverage = c(0.9, 0.95, 0.99), q = 0)

Arguments

Details

A target at or below a seed's final coverage is read off that seed's spatial accumulation curve by linear interpolation in coverage. Above it, every seed has accumulated the same pooled sample, so the value is the coverage-based extrapolation of the pooled abundances (Chao et al. 2014): the sample size reaching the target coverage is found from the Chao-Jost coverage deficit, and the Hill number at that size approaches the Chao1 (q = 0), Chao et al. (2013) Shannon (q = 1) or unbiased Simpson (q = 2) asymptote; target_coverage = 1 gives the asymptote itself. Extrapolation is defined for integer q and for the individual-based coverage (coverage = "chao").

Value

An object of class spacc_coverage_ext containing:

References

Chao, A. & Jost, L. (2012). Coverage-based rarefaction and extrapolation: standardizing samples by completeness rather than size. Ecology, 93, 2533-2547.

Chao, A., Gotelli, N.J., Hsieh, T.C., et al. (2014). Rarefaction and extrapolation with Hill numbers: a framework for sampling and estimation in species diversity studies. Ecological Monographs, 84, 45-67.

See Also

spaccCoverage(), interpolateCoverage()

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rpois(50 * 30, 2), nrow = 50)
cov <- spaccCoverage(species, coords)
ext <- extrapolateCoverage(cov, target_coverage = c(0.95, 0.99))
print(ext)

Gamma Diversity (Regional)

Description

Compute Hill numbers for the pooled community across all sites.

Usage

gammaDiversity(x, q = c(0, 1, 2))

Arguments

Details

Gamma diversity represents regional (total) diversity across all sites. It is computed by pooling abundances across all sites and calculating Hill numbers on the combined community.

Value

A named numeric vector with gamma diversity for each q.

References

Jost, L. (2007). Partitioning diversity into independent alpha and beta components. Ecology, 88, 2427-2439.

See Also

alphaDiversity() for local diversity, diversityPartition() for full alpha-beta-gamma decomposition

Examples

species <- matrix(rpois(50 * 30, 2), nrow = 50)
gammaDiversity(species, q = c(0, 1, 2))

Improved Chao1 (iChao1) Richness Estimator

Description

Estimate total species richness using the improved Chao1 estimator (Chiu et al. 2014). Uses singletons through quadrupletons (f1–f4) to reduce bias for small samples.

Usage

iChao1(x)

Arguments

Details

The improved Chao1 estimator adds a correction term using f3 and f4:

S_{iChao1} = S_{Chao1} + \frac{f_3}{4 f_4} \max\left(f_1 - \frac{f_2 f_3}{2 f_4}, 0\right)

When f_4 = 0, the estimator collapses to Chao1.

Value

An object of class spacc_estimate with components:

estimator

Name of the estimator ("iChao1")

estimate

Estimated total richness

se

Standard error of the estimate

lower

Lower 95 percent confidence bound

upper

Upper 95 percent confidence bound

S_obs

Observed species richness

details

List with f1, f2, f3, f4

References

Chiu, C.H., Wang, Y.T., Walther, B.A. & Chao, A. (2014). An improved nonparametric lower bound of species richness via a modified Good-Turing frequency formula. Biometrics, 70, 671-682.

See Also

chao1() for the standard estimator, iChao2() for incidence-based version

Examples

species <- matrix(rpois(50 * 30, 2), nrow = 50)
iChao1(species)

Improved Chao2 (iChao2) Richness Estimator

Description

Estimate total species richness from incidence data using the improved Chao2 estimator (Chiu et al. 2014). Uses uniques through quadruplicates (Q1–Q4) to reduce bias for small samples.

Usage

iChao2(x)

Arguments

Details

The improved Chao2 estimator is the incidence-based analogue of iChao1:

S_{iChao2} = S_{Chao2} + \frac{Q_3}{4 Q_4} \max\left(Q_1 - \frac{Q_2 Q_3}{2 Q_4}, 0\right)

When Q_4 = 0, the estimator collapses to Chao2.

Value

An object of class spacc_estimate with components:

estimator

Name of the estimator ("iChao2")

estimate

Estimated total richness

se

Standard error of the estimate

lower

Lower 95 percent confidence bound

upper

Upper 95 percent confidence bound

S_obs

Observed species richness

details

List with Q1, Q2, Q3, Q4, n_sites

References

Chiu, C.H., Wang, Y.T., Walther, B.A. & Chao, A. (2014). An improved nonparametric lower bound of species richness via a modified Good-Turing frequency formula. Biometrics, 70, 671-682.

See Also

chao2() for the standard estimator, iChao1() for abundance-based version

Examples

species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
iChao2(species)

Interpolate Richness at Target Coverage Levels

Description

Estimate species richness at specified coverage levels by interpolation.

Usage

interpolateCoverage(x, target = c(0.9, 0.95, 0.99))

Arguments

Value

A data.frame with columns for each target coverage level, showing interpolated richness for each seed.

References

Chao, A. & Jost, L. (2012). Coverage-based rarefaction and extrapolation: standardizing samples by completeness rather than size. Ecology, 93, 2533-2547.

Jackknife Richness Estimator

Description

First- or second-order jackknife estimator of total species richness (Burnham & Overton 1978, 1979).

Usage

jackknife(x, order = 1L)

Arguments

Details

First-order jackknife:

S_{jack1} = S_{obs} + Q_1 \frac{n-1}{n}

Second-order jackknife:

S_{jack2} = S_{obs} + \frac{Q_1(2n-3)}{n} - \frac{Q_2(n-2)^2}{n(n-1)}

Value

An object of class spacc_estimate.

References

Burnham, K.P. & Overton, W.S. (1978). Estimation of the size of a closed population when capture probabilities vary among animals. Biometrika, 65, 625-633.

Burnham, K.P. & Overton, W.S. (1979). Robust estimation of population size when capture probabilities vary among animals. Ecology, 60, 927-936.

See Also

chao2(), bootstrap_richness()

Examples

species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
jackknife(species, order = 1)
jackknife(species, order = 2)

Exact (Mao Tau) Expected Accumulation

Description

Compute the expected species accumulation curve using sample-based rarefaction (Mao Tau estimator). This is analytically identical to the expected curve from random permutations.

Usage

mao_tau(x)

Arguments

Value

A data.frame with columns: sites, expected, sd, lower, upper

References

Colwell, R.K., Mao, C.X. & Chang, J. (2004). Interpolating, extrapolating, and comparing incidence-based species accumulation curves. Ecology, 85, 2717-2727.

Plot Spatial SAC

Description

Create a ggplot2 visualization of species accumulation curves. For grouped spacc objects, curves are overlaid with different colors.

Usage

## S3 method for class 'spacc'
plot(
  x,
  ci = TRUE,
  ci_level = 0.95,
  ci_alpha = 0.3,
  show_seeds = FALSE,
  saturation = FALSE,
  saturation_level = 0.9,
  facet = FALSE,
  ...
)

Arguments

Value

A ggplot2 object.

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
sac <- spacc(species, coords)
plot(sac)

Plot spacc_metrics

Description

Create visualizations of per-site accumulation metrics.

Usage

## S3 method for class 'spacc_metrics'
plot(
  x,
  metric = NULL,
  type = c("heatmap", "points", "histogram"),
  point_size = 3,
  palette = "viridis",
  ...
)

Arguments

Value

A ggplot2 object.

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
metrics <- spaccMetrics(species, coords, metrics = c("slope_10", "auc"))
plot(metrics, metric = "slope_10", type = "heatmap")

Predict from Empirical Accumulation Curve

Description

Interpolate the mean empirical accumulation curve at arbitrary site counts using linear interpolation. Unlike the predict method for spacc_fit objects, this does not use a fitted model; it interpolates the observed curve directly.

Usage

## S3 method for class 'spacc'
predict(object, n = NULL, ci = TRUE, ci_level = 0.95, warn = TRUE, ...)

Arguments

Value

A data frame (if ci = TRUE) or named numeric vector (if ci = FALSE). Out-of-range values return NA.

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
sac <- spacc(species, coords, n_seeds = 10, progress = FALSE)
predict(sac, n = c(10, 25, 50))
predict(sac, n = c(10, 25), ci = FALSE)

Predict richness at a given area

Description

Predict richness at a given area

Usage

## S3 method for class 'spacc_area'
predict(
  object,
  area = NULL,
  interval = c("none", "bootstrap"),
  level = NULL,
  ...
)

Arguments

Value

A numeric vector (point) or a data frame with area, fit, lower, upper.

Predict richness at a given sampling effort

Description

Predict richness at a given sampling effort

Usage

## S3 method for class 'spacc_fit'
predict(object, n = NULL, interval = c("none", "bootstrap"), level = NULL, ...)

Arguments

Value

A numeric vector (point) or a data frame with n, fit, lower, upper (interval).

Individual-Based Rarefaction

Description

Compute individual-based rarefaction curves for Hill numbers at any order q. This complements the sample-based accumulation in spacc().

Usage

rarefy(x, n_individuals = NULL, q = 0, n_boot = 100L)

Arguments

Details

The expected Hill number of a random subsample of n individuals drawn without replacement from the pooled sample (Chao et al. 2014). Each species' subsample count is hypergeometric, and

{}^qD(n) = \left(\sum_i E\left[(Y_i / n)^q\right]\right)^{1/(1-q)},

with the exponential of the expected Shannon entropy at q = 1. At q = 0 this is the Hurlbert (1971) formula. Sizes above the sample size are extrapolated towards the asymptotic estimators (Chao et al. 2014; integer q only).

Value

An object of class spacc_rare containing:

References

Hurlbert, S.H. (1971). The nonconcept of species diversity: a critique and alternative parameters. Ecology, 52, 577-586.

Chao, A., Gotelli, N.J., Hsieh, T.C., et al. (2014). Rarefaction and extrapolation with Hill numbers: a framework for sampling and estimation in species diversity studies. Ecological Monographs, 84, 45-67.

Examples

abundance_matrix <- matrix(rpois(50 * 30, 2), nrow = 50)
rare <- rarefy(abundance_matrix)
plot(rare)

# Shannon rarefaction
rare_q1 <- rarefy(abundance_matrix, q = 1)
plot(rare_q1)

Standardized Effect Size (SES) via Null Models

Description

Compute standardized effect sizes by comparing observed diversity metrics against a null distribution from randomized communities. Supports multiple null model algorithms and works with most spacc output classes.

Usage

ses(
  x,
  species,
  coords = NULL,
  metric = NULL,
  null_model = c("frequency", "richness", "both", "curveball", "torus", "spatial_swap",
    "shift_rotate"),
  n_perm = 999L,
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  seed = NULL
)

Arguments

Details

SES is computed as:

SES = \frac{observed - \bar{null}}{sd_{null}}

A two-tailed p-value is calculated as the proportion of null values at least as extreme as the observed value:

p = \frac{2 \cdot \min(r, n_{perm} + 1 - r)}{n_{perm} + 1}

where r is the rank of the observed value among null values.

Null model algorithms:

• "frequency": Tests whether species composition matters given observed species frequencies

• "richness": Tests whether species identity matters given observed site richness

• "both": Maintains both marginal totals; tests non-random species co-occurrence patterns

• "curveball": Efficient alternative to "both" with proven uniform sampling properties

• "shift_rotate": The permutation nulls destroy the spatial structure of the community, so a spatially structured curve differs from almost any of them. The shift-and-rotate null keeps that structure intact and asks instead whether the curve depends on where the fixed design sits within it: each replicate moves the focal points by a rigid motion (rotation about their centroid, shift within the extent of the sites), drops those that leave the convex hull of the sites, and recomputes the curve. It is defined only for a fixed design; randomly drawn focal points are already a sample over all positions.

For method = "knn" analyses the permutation nulls reuse the observed focal points, so observed and null curves share one design.

Value

An object of class spacc_ses containing:

References

Gotelli, N.J. (2000). Null model analysis of species co-occurrence patterns. Ecology, 81, 2606-2621.

Strona, G., Nappo, D., Boccacci, F., Fattorini, S. & San-Miguel-Ayanz, J. (2014). A fast and unbiased procedure to randomize ecological binary matrices with fixed row and column totals. Nature Communications, 5, 4114.

Ridder, G., Hardy, O.J. & Ovaskainen, O. (2024). Generating spatially realistic environmental null models with the shift-&-rotate approach helps evaluate false positives in species distribution modelling. Methods in Ecology and Evolution, 15, 2331-2342. doi:10.1111/2041-210X.14443

See Also

spaccHill(), spaccBeta(), spaccMetrics()

Examples

coords <- data.frame(x = runif(20), y = runif(20))
species <- matrix(rbinom(20 * 15, 1, 0.3), nrow = 20)

sac <- spacc(species, coords, n_seeds = 10)
result <- ses(sac, species, n_perm = 19)
print(result)

Sampling Effort Species-Area Relationship (SESARS)

Description

Model the joint effect of sampling effort and area on species richness. Corrects for unequal survey intensity across sites, common in atlas data and citizen science datasets.

Usage

sesars(object, effort, model = c("power", "additive"), ...)

Arguments

Details

Standard SARs assume complete sampling within each area unit. SESARS incorporates sampling effort (E) alongside area (A) to provide unbiased richness estimates across regions with unequal survey intensity.

Value

An object of class spacc_sesars containing:

References

Dennstadt, F., Horak, J. & Martin, M.D. (2019). Predictive sampling effort and species-area relationship models for estimating richness in fragmented landscapes. Diversity and Distributions, 26, 1112-1123.

See Also

extrapolate(), spacc()

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
sac <- spacc(species, coords)
effort <- rpois(50, 10)
ses <- sesars(sac, effort, model = "power")
print(ses)

Species-Fragmented Area Relationship (SFAR)

Description

Separate the effects of habitat loss (area reduction) from fragmentation (splitting into patches) on species richness. Extends the classic power-law SAR with an explicit fragmentation term.

Usage

sfar(object, patches, model = c("power", "log"), ...)

Arguments

Details

The SFAR (Hanski et al. 2013) extends the power-law SAR to quantify the additional effect of habitat fragmentation on species richness. The model S = c * A^z * n^(-f) adds a penalty term for fragmentation (n = number of fragments), where f > 0 indicates that fragmentation reduces richness beyond what area loss alone would predict.

Value

An object of class spacc_sfar containing:

References

Hanski, I., Zurita, G.A., Bellocq, M.I. & Rybicki, J. (2013). Species-fragmented area relationship. Proceedings of the National Academy of Sciences, 110, 12715-12720.

Rybicki, J. & Hanski, I. (2013). Species-area relationships and extinctions caused by habitat loss and fragmentation. Ecology Letters, 16, 27-38.

See Also

extrapolate(), sesars()

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
sac <- spacc(species, coords)
patches <- kmeans(coords, centers = 5)$cluster
sfar_result <- sfar(sac, patches)
print(sfar_result)

Spatial Species Accumulation Curves

Description

Compute species accumulation curves using various spatial sampling methods with C++ backend for performance.

Usage

spacc(
  x,
  coords,
  n_seeds = 50L,
  method = c("knn", "kncn", "nn_walk", "random", "gaussian", "cone", "collector"),
  distance = c("euclidean", "haversine"),
  backend = c("auto", "exact", "kdtree"),
  support = NULL,
  include_halo = TRUE,
  sigma = NULL,
  cone_width = pi/4,
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  groups = NULL,
  time = NULL,
  w_space = 1,
  w_time = 1,
  seed = NULL,
  focal_points = NULL,
  focal_domain = NULL,
  order = NULL
)

Arguments

Details

The canonical knn method follows Chiarucci et al. (2009): each curve uses one fixed focal point and accumulates sites by increasing distance from that point. The nn_walk method is a greedy traversal whose reference point moves to the most recently selected site. The two methods represent different spatial sampling models.

Every method is a rule for ordering the sites; the species data are then accumulated along those orders. The orders element of the result lets the same traversals be reused, for example through the order argument.

Value

When groups = NULL, an object of class spacc containing:

References

Arrhenius, O. (1921). Species and area. Journal of Ecology, 9, 95-99.

Scheiner, S.M. (2003). Six types of species-area curves. Global Ecology and Biogeography, 12, 441-447.

Chiarucci, A., Bacaro, G., Rocchini, D., Ricotta, C., Palmer, M.W. & Scheiner, S.M. (2009). Spatially constrained rarefaction: incorporating the autocorrelated structure of biological communities into sample-based rarefaction. Community Ecology, 10, 209-214.

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)

# Basic usage
sac <- spacc(species, coords)
plot(sac)

# Different methods
sac_knn <- spacc(species, coords, method = "knn")
sac_rand <- spacc(species, coords, method = "random")
comp <- compare(sac_knn, sac_rand)

Spatial Beta Diversity Accumulation

Description

Analyze how beta diversity changes as sites are accumulated spatially. Partitions beta diversity into turnover (species replacement) and nestedness (species loss) components following Baselga (2010).

Usage

spaccBeta(
  x,
  coords,
  traits = NULL,
  tree = NULL,
  n_seeds = 50L,
  method = c("knn", "nn_walk"),
  index = c("sorensen", "jaccard"),
  distance = c("euclidean", "haversine"),
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  seed = NULL,
  map = FALSE,
  focal_points = NULL,
  focal_domain = NULL,
  pool = NULL
)

Arguments

Details

At each step of spatial accumulation, beta diversity is calculated between the accumulated species pool and the newly added site. This reveals how species composition changes as you expand spatially.

Interpretation:

• High turnover: New sites contribute different species (replacement)

• High nestedness: New sites contribute subsets of existing species (loss)

The sum of turnover and nestedness equals total beta diversity.

Supplying traits computes functional beta diversity, weighting the partition by the trait distinctiveness of the species exchanged between the accumulated pool and each new site. Supplying tree computes phylogenetic beta diversity (PhyloSor), weighting by shared branch length. The taxonomic, functional, and phylogenetic variants all return a spacc_beta object whose beta_type field records which was computed. map = TRUE is supported for the taxonomic variant only.

With pool, sites are aggregated into pooled assemblages (a species is present in a pool if it occurs in any of its sites) and the partition is repeated between pools: each site order becomes the order in which its pools are first reached, and the distance between two pools is the mean distance between their sites. Partitions at the two resolutions can disagree, for example when degraded habitats are homogenised at the pool level while individual sites still differ (Jones et al. 2025). summary(x, resolution = "pooled") and plot(x, resolution = "pooled") or "both" show the pooled partition.

Value

An object of class spacc_beta containing:

References

Baselga, A. (2010). Partitioning the turnover and nestedness components of beta diversity. Global Ecology and Biogeography, 19, 134-143.

Jones, F.A.M., Hardenbol, A.A., Hekkala, A.-M., et al. (2025). Partitioning beta diversity at two spatial resolutions reveals biotic homogenisation with habitat degradation. Diversity and Distributions, 31, e70080. doi:10.1111/ddi.70080

See Also

betapart::beta.pair() for pairwise beta diversity

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)

beta <- spaccBeta(species, coords, n_seeds = 30)
plot(beta)

# Compare Sorensen vs Jaccard
beta_jac <- spaccBeta(species, coords, index = "jaccard")

# Functional beta diversity (pass a species-by-traits matrix)
traits <- matrix(rnorm(30 * 3), nrow = 30)
rownames(traits) <- colnames(species) <- paste0("sp", 1:30)
beta_func <- spaccBeta(species, coords, traits = traits)

Coverage-Based Spatial Rarefaction

Description

Compute spatial accumulation curves with sample coverage tracking. Allows standardization by completeness (coverage) rather than sample size, following Chao & Jost (2012) or the sample-based estimator of Chiu (2023).

Usage

spaccCoverage(
  x,
  coords,
  n_seeds = 50L,
  method = c("knn", "nn_walk"),
  distance = c("euclidean", "haversine"),
  coverage = c("chao", "chiu"),
  n_units = NULL,
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  seed = NULL,
  map = FALSE,
  focal_points = NULL,
  focal_domain = NULL
)

Arguments

Details

Sample coverage estimates the proportion of the total community abundance represented by observed species. It provides a measure of sampling completeness that is independent of sample size.

The Chao-Jost (2012) estimator counts singletons (f1) and doubletons (f2) in the cumulative abundance vector. It assumes individuals are sampled independently, which may not hold for plot-based spatial data.

The Chiu (2023) estimator uses incidence frequency counts instead: Q1 (species in exactly 1 site), Q2 (species in exactly 2 sites), and G1 (total abundance of Q1 species). It gives near-unbiased coverage estimates when organisms are spatially aggregated across sampling units. It writes the abundance of undetected species as the Chao2 estimate of their number, \frac{t-1}{t} \frac{Q_1^2}{2 Q_2} after t sites, times the mean abundance of a unique species, G_1 / Q_1.

With n_units = T, the t sites are a fraction w = t / T of a finite set of units drawn without replacement, and the number of undetected species is the MoRE term of chao2(x, n_units = ) (Chiu 2023, Methods in Ecology and Evolution, Eq. 5),

\frac{t-1}{t} \frac{(1-w) Q_1^2}{Q_2 + w (Q_1 - Q_2)^+},

so richness and coverage of the same sample rest on one estimate of the undetected species. It vanishes at w = 1: the coverage reaches 1 once the accumulation has visited every unit. On simulated grids its error is lower than that of the with-replacement form from w = 0.5 on, and several-fold lower near a census; it is about equal at w = 0.3 and larger at w = 0.1 on a sparse grid (vignette("extrapolation")).

Coverage-based rarefaction allows fair comparison of diversity across communities with different abundances by standardizing to equal completeness rather than equal sample size.

Value

An object of class spacc_coverage containing:

References

Chao, A. & Jost, L. (2012). Coverage-based rarefaction and extrapolation: standardizing samples by completeness rather than size. Ecology, 93, 2533-2547.

Chiu, C.-H. (2023). A species richness estimator for sample-based incidence data sampled without replacement. Methods in Ecology and Evolution, 14, 2482-2493. doi:10.1111/2041-210X.14146

Chiu, C.-H. (2023). A sample-based estimator for sample coverage. Ecology, 104, e4099.

See Also

iNEXT::iNEXT() for coverage-based rarefaction without spatial structure

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rpois(50 * 30, 2), nrow = 50)

cov <- spaccCoverage(species, coords)
plot(cov)

# Sample-based coverage (recommended for spatial data)
cov_chiu <- spaccCoverage(species, coords, coverage = "chiu")

# The 50 sites are 80% of a 62-cell grid
cov_grid <- spaccCoverage(species, coords, coverage = "chiu", n_units = 62)

# Interpolate richness at 90% and 95% coverage
interp <- interpolateCoverage(cov, target = c(0.90, 0.95))

Spatial Accumulation of a Custom Diversity Metric

Description

Accumulate any user-supplied diversity index along a spatial ordering of sites. At each accumulation step the cumulative community is passed to fun, which returns a single number. This is the general escape hatch behind the built-in metric functions: use it for indices that spacc does not implement directly.

Usage

spaccDiversity(
  x,
  coords,
  fun,
  ...,
  method = c("knn", "kncn", "nn_walk", "random", "collector"),
  incidence = FALSE,
  n_seeds = 50L,
  distance = c("euclidean", "haversine"),
  progress = TRUE,
  seed = NULL,
  focal_points = NULL,
  focal_domain = NULL
)

Arguments

Details

The site ordering reuses the same spatial traversals as the built-in methods (fixed-focus, nearest-neighbour walk, nearest-centroid, random, or data order), then evaluates fun on the accumulating community. Because the index is an arbitrary R function, this trades the speed of the compiled metrics for full flexibility.

Value

An object of class spacc_diversity that inherits from spacc, so the standard summary(), plot(), as.data.frame() and predict() methods apply. curves is an ⁠n_seeds x n_sites⁠ matrix of the metric along the accumulation.

See Also

spaccHill(), spaccPhylo(), spaccFunc() for built-in metrics.

Examples

coords <- data.frame(x = runif(40), y = runif(40))
species <- matrix(rpois(40 * 20, 2), nrow = 40)

# Shannon entropy along the accumulation
shannon <- function(comm) {
  p <- comm[comm > 0] / sum(comm)
  -sum(p * log(p))
}
div <- spaccDiversity(species, coords, shannon, n_seeds = 20)
plot(div)

Spatial Endemism Accumulation

Description

Compute the number of endemic species (species found only within the accumulated area) as sites are added spatially. Complements the standard SAC by tracking species unique to each spatial extent.

Usage

spaccEndemism(
  x,
  coords,
  n_seeds = 50L,
  method = c("knn", "nn_walk"),
  distance = c("euclidean", "haversine"),
  map = FALSE,
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  seed = NULL,
  focal_points = NULL,
  focal_domain = NULL
)

Arguments

Details

At each accumulation step k, an endemic species is one that is present in the accumulated sites (1..k) but absent from all remaining unvisited sites (k+1..n). This tracks how many species are unique to the area sampled so far.

The endemism curve typically starts low (few endemics at small areas), increases as the region grows, and eventually equals total richness when all sites are included.

Value

An object of class spacc_endemism containing:

References

Kier, G., Kreft, H., Lee, T.M., et al. (2009). A global assessment of endemism and species richness across island and mainland regions. Proceedings of the National Academy of Sciences, 106, 9322-9327.

May, F., Gerstner, K., McGlinn, D.J., et al. (2018). mobsim: an R package for the simulation and measurement of biodiversity across spatial scales. Methods in Ecology and Evolution, 9, 1401-1408.

See Also

spacc(), spaccHill()

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)

end <- spaccEndemism(species, coords, n_seeds = 30)
plot(end)

Spatial Functional Diversity Accumulation

Description

Compute spatial accumulation of functional diversity metrics based on traits.

Usage

spaccFunc(
  x,
  coords,
  traits,
  metric = c("fdis", "fric"),
  n_seeds = 50L,
  method = c("knn", "nn_walk"),
  distance = c("euclidean", "haversine"),
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  seed = NULL,
  map = FALSE,
  focal_points = NULL,
  focal_domain = NULL,
  axis_length_k = NULL,
  itv = FALSE
)

Arguments

Details

Functional diversity metrics quantify trait space occupation:

• FDis (Functional Dispersion): Abundance-weighted mean distance from the community centroid in trait space. Captures functional divergence.

• FRic (Functional Richness): Volume of trait space occupied (convex hull). Requires more species than traits to compute.

• K (aggregated functional diversity): the geometric mean of four facets between 0 and 1, each equal to 1 when species are spread uniformly along the trait axis with equal abundances (Wojcik et al. 2024). Per trait, with the n present species sorted by trait value and relative abundances b_i: functional richness is the trait range over the axis length; biomass evenness compares the b_i with 1/n; trait evenness compares the gaps between neighbouring trait values with an even spacing; dispersion is the abundance-weighted trait variance over that of a uniform distribution on the same range, capped at 1. Several traits combine as the geometric mean over traits, and each reported facet is its geometric mean over traits. Trait evenness needs three species: with two, K is the cube root of the other three facets and k_trait_evenness is NA. With fewer than two species K is 0. With presence/absence data biomass evenness is 1.

Intraspecific trait variability (itv = TRUE; Palacio et al. 2025). Each trait row is a unit of the community in place of a species: at each step the metrics are evaluated on the rows describing the accumulated sites, each weighted by its species' abundance at its site, shared equally among the rows describing that species there. Rows without a site (no site column, or NA) stand for the species at every site where it has no site-specific rows, so one row per species reproduces the species-mean curves. A widespread species whose traits change along a gradient then adds trait range as the accumulation reaches new parts of the gradient, which its mean trait value hides.

Value

An object of class spacc_func containing:

References

Laliberté, E. & Legendre, P. (2010). A distance-based framework for measuring functional diversity from multiple traits. Ecology, 91, 299-305.

Palacio, F.X., Ottaviani, G., Mammola, S., Graco-Roza, C., de Bello, F. & Carmona, C.P. (2025). Integrating intraspecific trait variability in functional diversity: an overview of methods and a guide for ecologists. Ecological Monographs, 95, e70024. doi:10.1002/ecm.70024

Wojcik, L.A., Gaedke, U., van Velzen, E. & Klauschies, T. (2024). Measuring overall functional diversity by aggregating its multiple facets: functional richness, biomass evenness, trait evenness and dispersion. Methods in Ecology and Evolution, 16, 215-227. doi:10.1111/2041-210X.14470

See Also

FD::dbFD() for comprehensive functional diversity analysis

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rpois(50 * 20, 2), nrow = 50)

# Trait matrix (species x traits)
traits <- matrix(rnorm(20 * 3), nrow = 20)
rownames(traits) <- paste0("sp", 1:20)
colnames(species) <- rownames(traits)

func <- spaccFunc(species, coords, traits, metric = c("fdis", "fric"))
plot(func)

Spatial Accumulation with Hill Numbers

Description

Compute spatial species accumulation curves using Hill numbers (effective number of species) instead of raw richness. Hill numbers unify diversity measures: q=0 is richness, q=1 is exponential Shannon, q=2 is inverse Simpson.

Usage

spaccHill(
  x,
  coords,
  q = c(0, 1, 2),
  n_seeds = 50L,
  method = c("knn", "nn_walk"),
  distance = c("euclidean", "haversine"),
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  seed = NULL,
  map = FALSE,
  focal_points = NULL,
  focal_domain = NULL
)

Arguments

Details

Hill numbers (Chao et al. 2014) provide a unified framework for diversity measurement. Unlike raw richness (q=0), higher-order Hill numbers (q=1, q=2) down-weight rare species, providing different perspectives on diversity.

The spatial accumulation of Hill numbers can reveal scale-dependent diversity patterns missed by richness alone.

Value

An object of class spacc_hill containing:

References

Chao, A., Gotelli, N.J., Hsieh, T.C., Sander, E.L., Ma, K.H., Colwell, R.K. & Ellison, A.M. (2014). Rarefaction and extrapolation with Hill numbers: a framework for sampling and estimation in species diversity studies. Ecological Monographs, 84, 45-67.

See Also

spacc() for richness-only accumulation, iNEXT::iNEXT() for non-spatial Hill number rarefaction

Examples

# Compare diversity at different orders
coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rpois(50 * 30, 2), nrow = 50)

hill <- spaccHill(species, coords, q = c(0, 1, 2))
plot(hill)

# Extract summary at final site
summary(hill)

Spatial Hill Number Beta Diversity

Description

Compute multisite beta diversity as gamma/alpha decomposition of Hill numbers along the spatial accumulation curve (Jost 2007 framework).

Usage

spaccHillBeta(
  x,
  coords,
  q = c(0, 1, 2),
  n_seeds = 50L,
  distance = c("euclidean", "haversine"),
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  seed = NULL,
  method = c("knn", "nn_walk"),
  focal_points = NULL,
  focal_domain = NULL,
  standardize = c("none", "coverage"),
  coverage = NULL
)

Arguments

Details

At each accumulation step k, the function computes:

• Gamma: Hill number of the pooled community (all k sites combined)

• Alpha: Generalized mean of per-site Hill numbers (Jost's power mean)

• Beta: gamma / alpha (effective number of distinct communities)

Beta = 1 means all sites are identical; beta = k means all sites are completely different. This provides the Hill-number analogue of the Baselga-based spaccBeta().

Value

An object of class spacc_hill_beta containing:

Coverage standardisation

Observed beta depends on sampling completeness: gamma and alpha approach their asymptotes at different rates, so the ratio moves along the curve even where composition does not change. With standardize = "coverage", gamma is the Hill number of the pooled sample and alpha the Hill number of the joint sample (every species-by-site cell treated as a unit) divided by k (Chiu et al. 2014), both rarefied or extrapolated to the same sample coverage C with the estimators of Chao et al. (2014); beta is their ratio (Chao et al. 2023). By default C is, at each step, the lower of the two observed coverages, so both values come from the observed samples or their rarefaction. A fixed coverage makes the curve comparable across steps and datasets, at the cost of extrapolation where a step's sample is less complete; for q = 0, Chao et al. (2023) advise against extrapolating beyond the coverage of twice the sample size.

The standardised beta is transformed to four dissimilarities between 0 and 1 (Chao et al. 2023): 1-CqN (Sorensen-type, the classic Sorensen index at q = 0), 1-UqN (Jaccard-type), 1-VqN and 1-SqN. They are NA at the first step, where there is one site.

References

Jost, L. (2007). Partitioning diversity into independent alpha and beta components. Ecology, 88, 2427-2439.

Chiu, C.-H., Jost, L. & Chao, A. (2014). Phylogenetic beta diversity, similarity, and differentiation measures based on Hill numbers. Ecological Monographs, 84, 21-44.

Chao, A., Gotelli, N.J., Hsieh, T.C., et al. (2014). Rarefaction and extrapolation with Hill numbers: a framework for sampling and estimation in species diversity studies. Ecological Monographs, 84, 45-67.

Chao, A., Thorn, S., Chiu, C.-H., et al. (2023). Rarefaction and extrapolation with beta diversity under a framework of Hill numbers: the iNEXT.beta3D standardization. Ecological Monographs, 93, e1588. doi:10.1002/ecm.1588

See Also

spaccBeta() for P/A-based Baselga partitioning, spaccHill() for Hill accumulation without beta decomposition

Examples

coords <- data.frame(x = runif(40), y = runif(40))
species <- matrix(rpois(40 * 20, 2), nrow = 40)

hb <- spaccHillBeta(species, coords, n_seeds = 10, progress = FALSE)
plot(hb)

Spatial Hill Numbers at Standardized Coverage

Description

Compute spatial accumulation of Hill numbers alongside sample coverage, enabling standardized comparison at equal completeness levels (iNEXT-style analysis applied spatially).

Usage

spaccHillCoverage(
  x,
  coords,
  q = c(0, 1, 2),
  target_coverage = NULL,
  n_seeds = 50L,
  distance = c("euclidean", "haversine"),
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  seed = NULL,
  method = c("knn", "nn_walk"),
  focal_points = NULL,
  focal_domain = NULL
)

Arguments

Details

Combines spatial kNN accumulation with simultaneous Hill number and coverage computation. At each accumulation step, both the Chao-Jost sample coverage and Hill numbers for all requested orders are calculated.

When target_coverage is specified, Hill numbers are interpolated at those coverage levels using the existing interpolate_at_coverage() C++ function, enabling fair comparison between sites with different sampling completeness.

Value

An object of class spacc_hill_coverage containing:

References

Chao, A. & Jost, L. (2012). Coverage-based rarefaction and extrapolation: standardizing samples by completeness rather than size. Ecology, 93, 2533-2547.

Chao, A., Gotelli, N.J., Hsieh, T.C., Sander, E.L., Ma, K.H., Colwell, R.K. & Ellison, A.M. (2014). Rarefaction and extrapolation with Hill numbers. Ecological Monographs, 84, 45-67.

See Also

spaccHill() for Hill accumulation without coverage, spaccCoverage() for coverage accumulation without Hill numbers

Examples

coords <- data.frame(x = runif(40), y = runif(40))
species <- matrix(rpois(40 * 20, 2), nrow = 40)

hc <- spaccHillCoverage(species, coords, n_seeds = 10, progress = FALSE)
plot(hc)

# Standardize at 90% coverage
hc2 <- spaccHillCoverage(species, coords, target_coverage = 0.9,
                          n_seeds = 10, progress = FALSE)
summary(hc2)

Per-Site Accumulation Metrics

Description

Compute spatial accumulation metrics for each site as a starting point. Useful for identifying sites with high or low accumulation rates, visualizing spatial patterns in diversity, and understanding edge effects.

Usage

spaccMetrics(
  x,
  coords,
  metrics = c("slope_10", "half_richness", "auc"),
  method = c("knn", "kncn", "nn_walk", "random"),
  distance = c("euclidean", "haversine"),
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE
)

Arguments

Details

This function runs a spatial accumulation curve starting from each site individually, then extracts summary metrics from each curve. This allows you to identify:

• Sites in species-rich areas (high initial slope)

• Core vs edge sites (fast vs slow accumulation)

• Spatial patterns in community structure

The metrics can be plotted as a heatmap using plot(result, type = "heatmap"), which requires the ggplot2 package. For more sophisticated spatial visualization with study area boundaries, see the areaOfEffect package.

Value

An object of class spacc_metrics containing:

References

Soberon, J.M. & Llorente, J.B. (1993). The use of species accumulation functions for the prediction of species richness. Conservation Biology, 7, 480-488.

See Also

spacc() for standard accumulation curves

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
metrics <- spaccMetrics(species, coords,
                        metrics = c("slope_10", "auc"))
metrics$metrics$slope_10

Spatial Accumulation of Interaction Diversity

Description

Accumulate an interaction network over space: interaction richness and interaction Hill numbers, sample coverage, the richness of the two interacting levels (for example plants and pollinators) and connectance, all along the same site orders.

Usage

spaccNetwork(
  interactions,
  coords,
  q = c(0, 1, 2),
  n_seeds = 50L,
  method = c("knn", "nn_walk"),
  distance = c("euclidean", "haversine"),
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  seed = NULL,
  focal_points = NULL,
  focal_domain = NULL
)

Arguments

Details

Each unique lower-upper pair is a unit of diversity and its frequency its abundance, so the interaction Hill number of order q is the ordinary Hill number of the interaction frequencies of the accumulated sites, the effective number of equally frequent interactions, and the sample coverage is the Chao & Jost (2012) coverage of those frequencies (Chiu et al. 2023). The species curves of the two levels use each species' summed interaction frequency. Connectance is the realised share of the possible interactions between the species detected so far.

Value

An object of class spacc_network containing:

References

Chiu, C.-H., Chao, A., Vogel, S., Kriegel, P. & Thorn, S. (2023). Quantifying and estimating ecological network diversity based on incomplete sampling data. Philosophical Transactions of the Royal Society B, 378, 20220183. doi:10.1098/rstb.2022.0183

See Also

spaccHillCoverage() for the same quantities on species

Examples

set.seed(1)
coords <- data.frame(x = runif(30), y = runif(30))
interactions <- data.frame(
  site = sample(30, 400, replace = TRUE),
  lower = paste0("plant", sample(12, 400, replace = TRUE)),
  upper = paste0("bee", sample(20, 400, replace = TRUE))
)
net <- spaccNetwork(interactions, coords, n_seeds = 10)
net
plot(net)

Spatial Phylogenetic Diversity Accumulation

Description

Compute spatial accumulation of phylogenetic diversity metrics (MPD, MNTD, PD).

Usage

spaccPhylo(
  x,
  coords,
  tree,
  metric = c("mpd", "mntd"),
  rho = 0,
  n_seeds = 50L,
  method = c("knn", "nn_walk"),
  distance = c("euclidean", "haversine"),
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  seed = NULL,
  map = FALSE,
  focal_points = NULL,
  focal_domain = NULL
)

Arguments

Details

Phylogenetic diversity metrics incorporate evolutionary relationships:

• MPD (Mean Pairwise Distance): Average phylogenetic distance between all pairs of species. Sensitive to tree-wide patterns.

• MNTD (Mean Nearest Taxon Distance): Average distance to closest relative. Sensitive to terminal clustering.

• PD (Faith's Phylogenetic Diversity): Total branch length connecting species to the root. Requires full tree object.

• EvoHeritage: Features accumulate along each branch and are lost at rate rho, so a branch of length L passes on a fraction exp(-rho * L) of what entered it. A feature counts once if it survives to at least one present species. The value is in standardised units (the heritage at the end of a single branch of unit length), so rho = 0 gives PD and large rho tends to species richness. Scanning rho places a community between its richness and PD curves.

Tree-based metrics match the columns of x to tree$tip.label by name; the tree may contain species absent from x.

Value

An object of class spacc_phylo containing:

References

Faith, D.P. (1992). Conservation evaluation and phylogenetic diversity. Biological Conservation, 61, 1-10.

Rosindell, J., Manson, K., Gumbs, R., Pearse, W.D. & Steel, M. (2024). Phylogenetic biodiversity metrics should account for both accumulation and attrition of evolutionary heritage. Systematic Biology, 73, 158-182. doi:10.1093/sysbio/syad072

Webb, C.O. (2000). Exploring the phylogenetic structure of ecological communities: an example for rain forest trees. American Naturalist, 156, 145-155.

See Also

picante::mpd(), picante::mntd(), picante::pd()

Examples

if (requireNamespace("ape", quietly = TRUE)) {
  # Create random tree
  tree <- ape::rtree(30)

  coords <- data.frame(x = runif(50), y = runif(50))
  species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
  colnames(species) <- tree$tip.label

  phylo <- spaccPhylo(species, coords, tree, metric = c("mpd", "mntd"))
  plot(phylo)
}

Spatial Phylogenetic or Functional Redundancy Accumulation

Description

Track how much of the taxonomic diversity added along the spatial accumulation is new phylogenetic or functional diversity, as the ratio of the taxonomic Hill number to the phylogenetic or functional Hill number of the accumulated community.

Usage

spaccRedundancy(
  x,
  coords,
  tree = NULL,
  traits = NULL,
  q = c(0, 1, 2),
  n_seeds = 50L,
  method = c("knn", "nn_walk"),
  distance = c("euclidean", "haversine"),
  dist_method = c("euclidean", "gower"),
  normalize = TRUE,
  parallel = TRUE,
  n_cores = NULL,
  progress = TRUE,
  seed = NULL,
  focal_points = NULL,
  focal_domain = NULL
)

Arguments

Details

At each accumulation step the function computes, for every order q, the taxonomic Hill number ^qD of the accumulated community and either its phylogenetic Hill number ^qD(T) (Chao et al. 2010; effective number of equally abundant, equally distinct lineages, with T the abundance-weighted mean root-to-tip depth) or its similarity-sensitive functional Hill number ^qD^Z (Leinster & Cobbold 2012; with Z = 1 - d from trait distances). The ratio

R_q = {}^qD / {}^qD(T) \quad \textrm{or} \quad {}^qD / {}^qD^Z

is at least 1, equal to 1 when species are maximally distinct (a star tree, or every pair of species at the maximal trait distance), and grows as the species added are close relatives or functionally similar to those already present. A ratio rising along the curve says the area is adding species faster than lineages or trait combinations.

Value

An object of class spacc_redundancy containing:

References

Chao, A., Chiu, C.-H. & Jost, L. (2010). Phylogenetic diversity measures based on Hill numbers. Philosophical Transactions of the Royal Society B, 365, 3599-3609.

Leinster, T. & Cobbold, C.A. (2012). Measuring diversity: the importance of species similarity. Ecology, 93, 477-489.

Alberdi, A. & Gilbert, M.T.P. (2019). A guide to the application of Hill numbers to DNA-based diversity analyses. Molecular Ecology Resources, 19, 804-817. doi:10.1111/1755-0998.13014

See Also

spaccHill(), spaccPhylo(), spaccFunc(), diversityProfile()

Examples

if (requireNamespace("ape", quietly = TRUE)) {
  tree <- ape::rcoal(25)
  coords <- data.frame(x = runif(40), y = runif(40))
  species <- matrix(rpois(40 * 25, 1), nrow = 40,
                    dimnames = list(NULL, tree$tip.label))
  red <- spaccRedundancy(species, coords, tree = tree, n_seeds = 10)
  plot(red)
}

Wavefront Expansion Accumulation

Description

Accumulate species within an expanding radius from seed points. Models invasion spread from introduction points.

Usage

spaccWavefront(
  x,
  coords,
  n_seeds = 50L,
  r0 = 0,
  dr = NULL,
  n_steps = 50L,
  distance = c("euclidean", "haversine"),
  progress = TRUE,
  seed = NULL
)

Arguments

Value

An object of class spacc_wavefront containing:

References

Shigesada, N. & Kawasaki, K. (1997). Biological Invasions: Theory and Practice. Oxford University Press.

Examples

coords <- data.frame(x = runif(50), y = runif(50))
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)
wf <- spaccWavefront(species, coords, n_seeds = 20, n_steps = 50)
plot(wf)

Spatial Eigenvectors (PCNM/dbMEM)

Description

Compute spatial eigenvectors from site coordinates using Principal Coordinates of Neighbour Matrices (PCNM) or distance-based Moran's Eigenvector Maps (dbMEM).

Usage

spatialEigenvectors(
  coords,
  threshold = NULL,
  method = c("pcnm", "dbmem"),
  distance = c("euclidean", "haversine")
)

Arguments

Details

PCNM (Borcard & Legendre 2002) decomposes spatial structure into orthogonal eigenvectors representing patterns at different spatial scales. Large eigenvalues correspond to broad-scale patterns (positive spatial autocorrelation), while small eigenvalues represent fine-scale patterns.

The algorithm:

1. Compute pairwise distances

2. Truncate distances beyond threshold (set to 4 * threshold)

3. Double-centre the squared truncated distance matrix

4. Extract eigenvectors with positive eigenvalues

5. For "dbmem": additionally filter to Moran's I > 0

Value

An object of class spacc_mem containing:

References

Borcard, D. & Legendre, P. (2002). All-scale spatial analysis of ecological data by means of principal coordinates of neighbour matrices. Ecological Modelling, 153, 51-68.

Dray, S., Legendre, P. & Peres-Neto, P.R. (2006). Spatial modelling: a comprehensive framework for principal coordinate analysis of neighbour matrices (PCNM). Ecological Modelling, 196, 483-493.

See Also

spatialPartition() for variance partitioning with MEMs

Examples

coords <- data.frame(x = runif(30), y = runif(30))
mem <- spatialEigenvectors(coords)
print(mem)

Spatial Variance Partitioning with MEMs

Description

Partition diversity variation into spatial and non-spatial components using forward selection of Moran's Eigenvector Maps.

Usage

spatialPartition(x, mem, metric = NULL, forward = TRUE, alpha = 0.05)

Arguments

Details

Forward selection of MEMs proceeds by adding the MEM that most improves the model AIC at each step, stopping when no MEM improves AIC by more than 2 units or when p > alpha.

Value

An object of class spacc_mem_partition containing:

See Also

spatialEigenvectors() for computing MEMs

Examples

coords <- data.frame(x = runif(30), y = runif(30))
species <- matrix(rpois(30 * 15, 2), nrow = 30)

mem <- spatialEigenvectors(coords)
alpha <- alphaDiversity(species, q = 0)
part <- spatialPartition(alpha$q0, mem)
print(part)

Spatially-Constrained Rarefaction

Description

Compute expected species richness accounting for spatial autocorrelation (Chiarucci et al. 2009). Uses distance-weighted sampling probabilities.

Usage

spatialRarefaction(x, coords, n_perm = 100, bandwidth = NULL)

Arguments

Value

A data.frame with columns: sites, mean, sd, lower, upper

References

Chiarucci, A., Bacaro, G., Rocchini, D. & Fattorini, L. (2009). Discovering and rediscovering the sample-based rarefaction formula in the ecological literature. Community Ecology, 10, 195-199.

Spatial Subsampling

Description

Reduce spatial autocorrelation by subsampling sites using various methods.

Usage

subsample(
  coords,
  n = NULL,
  method = c("grid", "random", "thinning"),
  cell_size = NULL,
  min_dist = NULL,
  seed = NULL
)

Arguments

Details

Methods:

• "grid": Overlay a grid and retain one random site per cell.

• "random": Simple random subsample of n sites.

• "thinning": Iteratively remove sites until minimum distance is achieved.

Value

Integer vector of row indices to retain.

References

Aiello-Lammens, M.E., Boria, R.A., Radosavljevic, A., et al. (2015). spThin: an R package for spatial thinning of species occurrence records for use in ecological niche models. Ecography, 38, 541-545.

Lennon, J.J., Koleff, P., Greenwood, J.J.D. & Gaston, K.J. (2004). Contribution of rarity and commonness to patterns of species richness. Ecology Letters, 7, 81-87.

Examples

coords <- data.frame(x = runif(50) * 100, y = runif(50) * 100)
species <- matrix(rbinom(50 * 30, 1, 0.3), nrow = 50)

keep <- subsample(coords, method = "grid", cell_size = 20)
sac <- spacc(species[keep, ], coords[keep, ])

Zeta Diversity

Description

Compute zeta diversity — the mean number of species shared across k sites — for increasing orders of k. The zeta decline curve reveals community assembly processes: exponential decline suggests stochastic assembly, while power-law decline indicates niche-based assembly.

Usage

zetaDiversity(
  x,
  coords,
  orders = 1:10,
  n_samples = 100L,
  method = c("knn", "random"),
  distance = c("euclidean", "haversine"),
  seed = NULL,
  progress = TRUE
)

Arguments

Details

Zeta diversity of order k (\zeta_k) is the mean number of species shared across k sites. Key properties:

• \zeta_1 = mean species richness per site

• \zeta_2 = mean number of species shared by any two sites

• \zeta_k decreases monotonically with k

The zeta decline ratio (\zeta_k / \zeta_{k-1}) is diagnostic:

• Constant ratio: exponential decline (stochastic assembly)

• Increasing ratio: power-law decline (deterministic/niche-based assembly)

The knn method selects spatially nearest k sites from each focal site, which is ecologically meaningful for testing spatial turnover. The random method samples random k-site combinations, providing a null expectation.

Value

An object of class spacc_zeta containing:

References

Hui, C. & McGeoch, M.A. (2014). Zeta diversity as a concept and metric that unifies incidence-based biodiversity patterns. The American Naturalist, 184, 684-694.

Latombe, G., McGeoch, M.A., Nipperess, D.A. & Hui, C. (2018). zetadiv: an R package for computing compositional change across multiple sites, assemblages or cases. bioRxiv, 324897.

See Also

spaccBeta() for pairwise beta diversity, distanceDecay() for distance-decay relationships

Examples

coords <- data.frame(x = runif(30), y = runif(30))
species <- matrix(rbinom(30 * 20, 1, 0.3), nrow = 30)
zeta <- zetaDiversity(species, coords, orders = 1:5, n_samples = 50)
plot(zeta)