Nothing
#' Synoptic tables and calculation of group-wise frequencies, fidelity and
#' differential species character (wide species-sample matrix)
#'
#' @keywords internal
#' @noRd
syntable_wide <- function(matrix, groups, abund = "percentage", type = "percfreq", digits = NULL,
phi_method = "default",
phi_standard = c("none", "target", "all"),
phi_target_size = NULL, phi_alpha = NULL) {
phi_standard <- match.arg(phi_standard)
# Rounding: 3 for phi when 'digits' not supplied; 0 otherwise
if (is.null(digits)) {
digits_nonphi <- 0L
digits_phi <- 3L
} else {
d <- as.integer(digits)
digits_nonphi <- d
digits_phi <- d
}
# --- Input & cleaning -------------------------------------------------------
X <- as.matrix(matrix)
if (ncol(X) < 1 || nrow(X) < 1)
stop("Input 'matrix' must have >= 1 row (samples) and >= 1 column (species).")
if (is.null(colnames(X)))
stop("Column names (species) must be set.")
# Align 'groups' to sample order (row names optional)
if (any(is.na(groups))) stop("NA values in 'groups' are not allowed.")
if (!is.null(names(groups))) {
if (!is.null(rownames(X))) {
# align by sample names
if (!all(rownames(X) %in% names(groups)))
stop("Named 'groups' must include all sample (row) names of 'matrix'.")
groups <- groups[rownames(X)]
} else {
# no row names: fall back to position
if (length(groups) != nrow(X))
stop("Length of 'groups' must equal number of rows (samples).")
groups <- unname(groups)
message("Row names are missing; ignoring names(groups) and aligning by row order.")
}
} else {
if (length(groups) != nrow(X))
stop("Length of 'groups' must equal the number of rows (samples).")
}
grp_factor <- factor(groups)
grp_levels <- levels(grp_factor)
if (length(grp_levels) < 1) stop("No grp_levels in 'groups'.")
# Handle empty strings and NAs; robust PA conversion
Xc <- matrix(as.character(X), nrow = nrow(X), ncol = ncol(X),
dimnames = dimnames(X))
if (any(Xc == "", na.rm = TRUE)) Xc[Xc == ""] <- NA
Xnum <- suppressWarnings(matrix(as.numeric(Xc), nrow = nrow(X), ncol = ncol(X),
dimnames = dimnames(X)))
nonnum <- is.na(Xnum) & !is.na(Xc)
if (any(nonnum)) {
pres <- matrix(FALSE, nrow = nrow(X), ncol = ncol(X), dimnames = dimnames(X))
pres[!is.na(Xnum)] <- Xnum[!is.na(Xnum)] != 0
pres[nonnum] <- Xc[nonnum] != "0"
X <- matrix(0, nrow = nrow(X), ncol = ncol(X), dimnames = dimnames(X))
X[pres] <- 1
abund <- "pa"
message("Non-numeric abundance values detected -> using presence/absence scale.")
} else {
Xnum[is.na(Xnum)] <- 0
X <- Xnum
}
if (!abund %in% c("percentage", "pa"))
stop("Argument 'abund' must be either 'percentage' (cover) or 'pa' (presence/absence).")
# Basic sizes and group indicator
N <- nrow(X); p <- ncol(X)
Z <- stats::model.matrix(~ grp_factor - 1)
colnames(Z) <- grp_levels
b_vec <- colSums(Z) # group sizes
# Helper: finalize p x G numeric matrix -> data.frame with species rows
.finalize <- function(M_pg) {
df <- as.data.frame(M_pg)
rownames(df) <- colnames(X) # species
colnames(df) <- grp_levels
df
}
# --- Calculations by 'type' -------------------------------------------------
if (type == "totalfreq") {
P <- X > 0
counts_Gp <- crossprod(Z, P) # G x p
results <- list("syntable" = .finalize(t(counts_Gp)),
"samplesize" = b_vec)
} else if (type == "percfreq") {
P <- X > 0
counts_Gp <- crossprod(Z, P)
percfreq <- t(counts_Gp) # p x G
percfreq <- sweep(percfreq, 2, b_vec, "/", check.margin = FALSE) * 100
percfreq[, b_vec == 0] <- NA # avoid div/0 artefacts
percfreq <- round(percfreq, digits = digits_nonphi)
results <- list("syntable" = .finalize(percfreq),
"samplesize" = b_vec)
} else if (type == "mean" && abund == "percentage") {
sums_Gp <- crossprod(Z, X)
means <- t(sweep(sums_Gp, 1, b_vec, "/", check.margin = FALSE))
means[, b_vec == 0] <- NA
means <- round(means, digits = digits_nonphi)
results <- list("syntable" = .finalize(means),
"samplesize" = b_vec)
} else if (type == "mean" && abund == "pa") {
stop("Cannot calculate mean abundance in groups with presence/absence values.")
} else if (type == "median" && abund == "percentage") {
med <- matrix(NA, nrow = p, ncol = length(grp_levels))
for (j in seq_along(grp_levels)) {
idx <- which(grp_factor == grp_levels[j])
med[, j] <- if (length(idx)) apply(X[idx, , drop = FALSE], 2, stats::median) else NA
}
med <- round(med, digits = digits_nonphi)
results <- list("syntable" = .finalize(med),
"samplesize" = b_vec)
} else if (type == "median" && abund == "pa") {
stop("Cannot calculate median abundance in groups with presence/absence values.")
} else if (type == "diffspec") {
# needs >= 3 grp_levels and no empty grp_levels
if (length(grp_levels) < 3)
stop("'diffspec' requires at least 3 grp_levels.")
if (any(b_vec == 0))
stop("At least one group has zero samples; cannot run 'diffspec' (division by zero).")
P <- X > 0
counts_Gp <- crossprod(Z, P) # G x p
syn <- t(counts_Gp) # p x G
syn <- sweep(syn, 2, b_vec, "/", check.margin = FALSE) * 100
syn <- round(syn, digits = digits_nonphi)
pb <- utils::txtProgressBar(min = 1, max = 23, style = 3)
small <- matrix(NA, nrow = nrow(syn), ncol = ncol(syn), dimnames = dimnames(syn))
for (i in seq_len(nrow(syn))) small[i, ] <- sort(as.numeric(syn[i, ]))
utils::setTxtProgressBar(pb, 3)
Avsmall <- matrix(NA, nrow = nrow(syn), ncol = ncol(syn))
Avsmall[, 1] <- small[, 1]
for (i in seq_len(nrow(syn))) for (k in 2:ncol(syn)) {
Avsmall[i, k] <- mean(as.numeric(small[i, 1:k]))
}
utils::setTxtProgressBar(pb, 4)
CondAvsmall <- matrix("", nrow = nrow(syn), ncol = ncol(syn))
for (k in seq_len(nrow(syn))) {
CondAvsmall[k, 2] <- if (small[k, 2] >= 2 * Avsmall[k, 1] + 20) "1" else "0"
}
for (i in 3:ncol(syn)) for (k in seq_len(nrow(syn))) {
CondAvsmall[k, i] <- if (CondAvsmall[k, i - 1] == "1") "1"
else if (small[k, i] >= 2 * Avsmall[k, i - 1] + 20) "1" else "0"
}
utils::setTxtProgressBar(pb, 5)
Avsmall2 <- matrix(NA, nrow = nrow(syn), ncol = ncol(syn))
Avsmall2[, 1] <- Avsmall[, 1]
for (i in 2:ncol(syn)) for (k in seq_len(nrow(syn))) {
Avsmall2[k, i] <- if (CondAvsmall[k, i] == "1") Avsmall2[k, i - 1] else Avsmall[k, i]
}
utils::setTxtProgressBar(pb, 6)
large <- matrix(NA, nrow = nrow(syn), ncol = ncol(syn), dimnames = dimnames(syn))
for (i in seq_len(nrow(syn))) large[i, ] <- sort(as.numeric(syn[i, ]), decreasing = TRUE)
utils::setTxtProgressBar(pb, 7)
Avlarge <- matrix(NA, nrow = nrow(syn), ncol = ncol(syn))
Avlarge[, 1] <- large[, 1]
for (i in seq_len(nrow(syn))) for (k in 2:ncol(syn)) {
Avlarge[i, k] <- mean(as.numeric(large[i, 1:k]))
}
utils::setTxtProgressBar(pb, 8)
CondAvlarge <- matrix("", nrow = nrow(syn), ncol = ncol(syn))
for (k in seq_len(nrow(syn))) {
CondAvlarge[k, 2] <- if (large[k, 2] <= (Avlarge[k, 1] / 2) - 10) "1" else "0"
}
for (i in 3:ncol(syn)) for (k in seq_len(nrow(syn))) {
CondAvlarge[k, i] <- if (CondAvlarge[k, i - 1] == "1") "1"
else if (large[k, i] <= (Avlarge[k, i - 1] / 2) - 10) "1" else "0"
}
utils::setTxtProgressBar(pb, 9)
Avlarge2 <- matrix(NA, nrow = nrow(syn), ncol = ncol(syn))
Avlarge2[, 1] <- Avlarge[, 1]
for (i in 2:ncol(syn)) for (k in seq_len(nrow(syn))) {
Avlarge2[k, i] <- if (CondAvlarge[k, i] == "1") Avlarge2[k, i - 1] else Avlarge[k, i]
}
utils::setTxtProgressBar(pb, 10)
CondAvsmall2 <- matrix("", nrow = nrow(syn), ncol = ncol(syn))
for (i in 2:ncol(syn)) for (k in seq_len(nrow(syn))) {
CondAvsmall2[k, i] <- if (small[k, i] > (Avlarge2[k, i - 1] / 2) - 10) "1" else "0"
}
utils::setTxtProgressBar(pb, 11)
Avsmall3 <- matrix(NA, nrow = nrow(syn), ncol = ncol(syn))
Avsmall3[, 1] <- Avsmall2[, 1]
for (i in 2:ncol(syn)) for (k in seq_len(nrow(syn))) {
Avsmall3[k, i] <- if (CondAvsmall2[k, i] == "1") Avsmall3[k, i - 1] else Avsmall2[k, i]
}
utils::setTxtProgressBar(pb, 12)
CondAvlarge2 <- matrix("", nrow = nrow(syn), ncol = ncol(syn))
for (i in 2:ncol(syn)) for (k in seq_len(nrow(syn))) {
CondAvlarge2[k, i] <- if (large[k, i] < 2 * Avsmall2[k, i - 1] + 20) "1" else "0"
}
utils::setTxtProgressBar(pb, 13)
Avlarge3 <- matrix(NA, nrow = nrow(syn), ncol = ncol(syn))
Avlarge3[, 1] <- Avlarge2[, 1]
for (i in 2:ncol(syn)) for (k in seq_len(nrow(syn))) {
Avlarge3[k, i] <- if (CondAvlarge2[k, i] == "1") Avlarge3[k, i - 1] else Avlarge2[k, i]
}
utils::setTxtProgressBar(pb, 14)
Small_Avsmall3 <- matrix("", nrow = nrow(syn), ncol = ncol(syn))
for (i in seq_len(ncol(syn))) for (k in seq_len(nrow(syn))) {
if (is.na(syn[k, 1])) {
Small_Avsmall3[k, i] <- ""
} else if (syn[k, i] >= 2 * max(Avsmall3[k, ], na.rm = TRUE) + 20) {
Small_Avsmall3[k, i] <- "p"
} else if (syn[k, i] < 2 * max(Avsmall3[k, ], na.rm = TRUE) + 20) {
Small_Avsmall3[k, i] <- "n"
} else stop("Calculation of differential species failed (Small_Avsmall3).")
}
utils::setTxtProgressBar(pb, 15)
Large_AvLarge3 <- matrix("", nrow = nrow(syn), ncol = ncol(syn))
for (i in seq_len(ncol(syn))) for (k in seq_len(nrow(syn))) {
if (is.na(syn[k, 1])) {
Large_AvLarge3[k, i] <- ""
} else if (syn[k, i] > min(Avlarge3[k, ], na.rm = TRUE) / 2 - 10) {
Large_AvLarge3[k, i] <- "p"
} else if (syn[k, i] <= min(Avlarge3[k, ], na.rm = TRUE) / 2 - 10) {
Large_AvLarge3[k, i] <- "n"
} else stop("Calculation of differential species failed (Large_AvLarge3).")
}
utils::setTxtProgressBar(pb, 16)
Sum_cond <- matrix("", nrow = nrow(syn), ncol = ncol(syn))
for (i in seq_len(ncol(syn))) for (k in seq_len(nrow(syn))) {
Sum_cond[k, i] <- paste0(Small_Avsmall3[k, i], Large_AvLarge3[k, i])
}
utils::setTxtProgressBar(pb, 17)
Results <- matrix("-", nrow = nrow(syn), ncol = ncol(syn), dimnames = dimnames(syn))
for (i in seq_len(ncol(syn))) for (k in seq_len(nrow(syn))) {
Results[k, i] <- switch(Sum_cond[k, i],
"pp" = "p", "nn" = "n", "pn" = "pn", "-")
}
utils::setTxtProgressBar(pb, 18)
subresult <- matrix(0L, nrow = nrow(syn), ncol = ncol(syn))
subresult[Results %in% c("n", "p", "pn")] <- 1L
onlydiff <- Results[rowSums(subresult) > 0, , drop = FALSE]
others <- rownames(Results)[rowSums(subresult) == 0]
# Per-group lists
positivespec <- matrix(NA, nrow = length(grp_levels), ncol = 200)
negativespec <- matrix(NA, nrow = length(grp_levels), ncol = 200)
positivenegativespec <- matrix(NA, nrow = length(grp_levels), ncol = 200)
for (i in seq_along(grp_levels)) {
sp <- rownames(Results)
pos <- sp[Results[, i] == "p"]; neg <- sp[Results[, i] == "n"]; pn <- sp[Results[, i] == "pn"]
if (length(pos)) positivespec[i, seq_along(pos)] <- pos
if (length(neg)) negativespec[i, seq_along(neg)] <- neg
if (length(pn)) positivenegativespec[i, seq_along(pn)] <- pn
}
diffspeclist <- setNames(vector("list", length(grp_levels)), grp_levels)
for (i in seq_along(grp_levels)) {
pos <- positivespec[i, ][!is.na(positivespec[i, ])]
neg <- negativespec[i, ][!is.na(negativespec[i, ])]
pn <- positivenegativespec[i, ][!is.na(positivenegativespec[i, ])]
diffspeclist[[i]] <- list("positive diff" = if (length(pos)) pos else noquote("no positive diff species"),
"negative diff" = if (length(neg)) neg else noquote("no negative diff species"),
"positive/negative diff" = if (length(pn)) pn else noquote("no positive/negative diff species"))
}
utils::setTxtProgressBar(pb, 23); close(pb)
results <- list("syntable" = .finalize(Results),
"onlydiff" = as.data.frame(onlydiff),
"others" = others,
"samplesize" = b_vec,
"differentials" = diffspeclist)
} else if (type == "phi") {
# ----- Fidelity (binary) with group-size equalization Tichý & Chytrý (2006)
if (!phi_method %in% c("default","ochiai","uvalue"))
stop("Unknown phi_method. Use 'default', 'ochiai', or 'uvalue'.")
P <- X > 0
a_pg <- t(crossprod(Z, P)) # species x groups (presences per group)
n_vec <- rowSums(a_pg) # species totals
b_vec <- colSums(Z) # group sizes
Nv <- nrow(X)
G <- ncol(a_pg)
# one-tailed Fisher test on original counts
# If phi_alpha is provided, compute p = P[X >= a] with X ~ Hypergeom(N, b, n)
# using the original (non-equalized) 2x2 counts, and later zero-out cells with p >= alpha.
do_fisher <- !is.null(phi_alpha)
if (do_fisher) {
alpha <- as.numeric(phi_alpha)
if (!is.finite(alpha) || alpha <= 0 || alpha >= 1)
stop("phi_alpha must be a number in (0, 1).")
# p-values per species (rows) × group (cols) from original 2x2 tables
p_mat <- matrix(1, nrow = nrow(a_pg), ncol = length(b_vec),
dimnames = list(colnames(X), grp_levels))
for (j in seq_along(b_vec)) {
# one-tailed: P[X >= a] for X ~ Hypergeom(Nv, m=b_vec[j], k=n_vec)
p_mat[, j] <- stats::phyper(q = a_pg[, j] - 1L,
m = b_vec[j],
n = Nv - b_vec[j],
k = n_vec,
lower.tail = FALSE)
}
}
# r_in = a / b (species × groups), Σ r_in across groups per species
Rin <- sweep(a_pg, 2, b_vec, "/")
Rin[!is.finite(Rin)] <- 0
sumR <- rowSums(Rin)
if (identical(phi_standard, "all")) {
# Target group gets s_t percent; other groups equal share of the remaining percent
s_pct <- if (is.null(phi_target_size)) 100 / G else as.numeric(phi_target_size)
if (!is.finite(s_pct) || s_pct <= 0 || s_pct >= 100)
stop("phi_target_size must be a percentage in (0, 100).")
if (G < 2) stop("phi_standard = 'all' requires at least 2 groups.")
s_t <- (s_pct / 100) * Nv # conceptual size of the target group
s_o <- ((100 - s_pct) / 100) * Nv / (G - 1) # size of each non-target group
N_eq <- s_t + (G - 1) * s_o # total conceptual plots (= Nv)
# a'_ij and n'_ij for each species i & target column j
Aeq_all <- Rin * s_t
# n'_ij = s_t * r_in_ij + s_o * (Σ r_in_i• - r_in_ij)
Neq_all <- matrix(s_o * sumR, nrow = nrow(Rin), ncol = G) + (s_t - s_o) * Rin
if (phi_method == "ochiai") {
den <- sqrt(Neq_all * s_t)
den[den == 0] <- Inf
phitab <- Aeq_all / den
} else {
num <- N_eq * Aeq_all - Neq_all * s_t
den <- sqrt(Neq_all * (N_eq - Neq_all) * s_t * (N_eq - s_t))
den[den == 0] <- Inf
phi_mat <- num / den
phitab <- if (phi_method == "uvalue") phi_mat * sqrt(max(N_eq - 1, 0)) else phi_mat
}
} else if (identical(phi_standard, "target")) {
# Equalize ONLY the target group to s_t; others left un-equalized (pooled outside)
s_pct <- if (is.null(phi_target_size)) 100 / G else as.numeric(phi_target_size)
if (!is.finite(s_pct) || s_pct <= 0 || s_pct >= 100)
stop("phi_target_size must be a percentage in (0, 100).")
s_t <- (s_pct / 100) * Nv
# r_out per target column: (n - a)/(N - b)
Rout <- sweep(a_pg, 2, b_vec, function(a, b) {
num <- (n_vec - a); den <- (Nv - b)
ifelse(den > 0, num / den, 0)
})
Aeq <- Rin * s_t # a'
Neq <- Aeq + Rout * (Nv - s_t) # n'
if (phi_method == "ochiai") {
den <- sqrt(Neq * s_t); den[den == 0] <- Inf
phitab <- Aeq / den
} else {
num <- Nv * Aeq - Neq * s_t
den <- sqrt(Neq * (Nv - Neq) * s_t * (Nv - s_t))
den[den == 0] <- Inf
phi_mat <- num / den
phitab <- if (phi_method == "uvalue") phi_mat * sqrt(max(Nv - 1, 0)) else phi_mat
}
} else { # phi_standard == "none"
if (phi_method == "ochiai") {
phitab <- a_pg / sqrt(n_vec %o% b_vec)
} else {
num <- Nv * a_pg - n_vec %o% b_vec
den <- sqrt((n_vec * (Nv - n_vec)) %o% (b_vec * (Nv - b_vec)))
den[!is.finite(den) | den == 0] <- Inf
phi_mat <- num / den
phitab <- if (phi_method == "uvalue") phi_mat * sqrt(max(Nv - 1, 0)) else phi_mat
}
}
# zero-out non-significant cells if Fisher test requested
if (do_fisher) {
phitab[p_mat >= alpha] <- 0
}
# Clean, name, round, package result
phitab[!is.finite(phitab)] <- 0
rownames(phitab) <- colnames(X)
colnames(phitab) <- grp_levels
phitab <- round(phitab, digits = digits_phi)
results <- list("syntable" = as.data.frame(phitab),
"samplesize" = b_vec)
} else {
stop("Cannot calculate synoptic table. Check 'abund' (percentage|pa) and 'type' (totalfreq|percfreq|mean|median|diffspec|phi).")
}
return(invisible(results))
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.