From 9ee40fe8e2c735dc10bbe78d8d58c7f41941ceec Mon Sep 17 00:00:00 2001 From: theHumanBorch Date: Fri, 26 Jun 2026 07:10:33 -0500 Subject: [PATCH 1/4] Input scaling, auto bandwidth, adaptive epsilon scale= (none/zscore/robust/arcsinh) learned at fit, reapplied at predict. affinityParams$alpha="auto" sets RBF/Laplace bandwidth by the median heuristic. epsilonQuantile= adapts suppression to the antibody spread. Result gains antibodies_unscaled; predict() now reads result$antibodies so test labels match training. Matured antibodies seed a Lloyd refinement coerced to K: agglomerate when the network kept >K prototypes, k-means++ split when = 1) + targetK <- as.integer(targetK) + } + if (!is.null(epsilonQuantile)) { + stopifnot("epsilonQuantile must be in (0, 1)" = + is.numeric(epsilonQuantile) && length(epsilonQuantile) == 1L && + epsilonQuantile > 0 && epsilonQuantile < 1) + } # --- Affinity/distance metric guard --- # Clonal selection matures antibodies to maximize `affinityFunc`, but @@ -141,6 +192,13 @@ AINet <- R6::R6Class( "random_uniform", "kmeans++")), consolidate = isTRUE(consolidate), consolidationSteps = as.integer(consolidationSteps), + scale = scale, + scaleCofactor = scaleCofactor, + targetK = targetK, + epsilonQuantile = epsilonQuantile, + coverageBoost = isTRUE(coverageBoost), + coverageQuantile = coverageQuantile, + scaling = NULL, # populated at fit() with learned stats verbose = verbose ) @@ -178,6 +236,19 @@ AINet <- R6::R6Class( task <- match.arg(task, c("clustering", "classification")) X <- as.matrix(X) + + # ================================ + # 0. Input scaling (C) + # ================================ + # Learn the transform on this training matrix and stash the stats so + # predict() applies the identical map to new data. epsilon, mutation + # scale and every distance are in feature units, so a fixed default + # only behaves sensibly once features share a comparable scale. + scaling <- private$.learn_scaling(X, self$config$scale, + self$config$scaleCofactor) + X <- private$.apply_scaling(X, scaling) + self$config$scaling <- scaling + n <- nrow(X) d <- ncol(X) cfg <- self$config @@ -223,7 +294,18 @@ AINet <- R6::R6Class( # Affinity/distance params (base values; iter_alpha may be modulated # by ClassSwitcher below). + # alpha = "auto" sets the RBF/Laplace bandwidth by the median heuristic: + # alpha = 1 / median(||x_i - x_j||^2) over a subsample of the (scaled) + # data, so the kernel resolves structure at the data's own length scale + # instead of the fixed alpha = 1 that is arbitrary once features are + # scaled. Only meaningful for distance-based kernels (gaussian/laplace). base_alpha <- cfg$affinityParams$alpha %||% 1 + if (is.character(base_alpha) && identical(base_alpha, "auto")) { + base_alpha <- private$.auto_bandwidth(X, cfg$affinityFunc) + if (cfg$verbose) { + cat(sprintf("Auto bandwidth: alpha = %.5g\n", base_alpha)) + } + } c_p <- cfg$affinityParams$c %||% 1 p_p <- cfg$affinityParams$p %||% 2 iter_alpha <- base_alpha @@ -478,8 +560,18 @@ AINet <- R6::R6Class( # (c) Network suppression [C++] # Removes near-duplicate antibodies within an epsilon-ball in distFunc. + # With epsilonQuantile set, the threshold tracks the current antibody + # spread (a low quantile of their pairwise distances) so suppression is + # scale-free and adapts as the repertoire contracts, rather than using a + # fixed feature-unit epsilon that means different things per dataset. + eps_iter <- cfg$epsilon + if (!is.null(cfg$epsilonQuantile)) { + eps_iter <- private$.adaptive_epsilon( + self$repertoire$cells, cfg$distFunc, cfg$epsilonQuantile, + p_p, Sigma_inv, fallback = cfg$epsilon) + } keep <- network_suppression_cpp( - self$repertoire$cells, cfg$distFunc, cfg$epsilon, + self$repertoire$cells, cfg$distFunc, eps_iter, p_p, Sigma_inv ) kept_idx <- which(keep) @@ -560,11 +652,48 @@ AINet <- R6::R6Class( m <- nrow(A_final) } + # ================================ + # 4a'. Coverage boost (D) + # ================================ + # Clonal expansion follows density, so rare populations can end the run + # with no nearby prototype. Find the worst-covered points (lowest max + # affinity to any surviving antibody) and drop k-means++ seeds among them + # so the final assignment / forced-K refinement has candidate prototypes + # for those regions. Clustering only (new seeds carry no class label). + if (task == "clustering" && isTRUE(cfg$coverageBoost) && nrow(A_final) >= 1L) { + cov_aff <- compute_affinity_matrix(X, A_final, cfg$affinityFunc, + iter_alpha, c_p, p_p) + coverage <- apply(cov_aff, 1, max) + thr <- stats::quantile(coverage, cfg$coverageQuantile, names = FALSE, + na.rm = TRUE) + under <- which(coverage <= thr) + if (length(under) >= 2L) { + n_new <- min(length(under), max(2L, round(0.5 * cfg$nAntibodies))) + new_seeds <- init_kmeanspp_cpp(X[under, , drop = FALSE], n_new) + A_final <- rbind(A_final, new_seeds) + self$repertoire$cells <- A_final + if (cfg$verbose) { + cat(sprintf("Coverage boost: +%d antibodies for %d under-covered points\n", + n_new, length(under))) + } + } + } + # ================================ # 4b. Final assignment [C++] # ================================ if (task == "clustering") { - if (cfg$consolidate && cfg$consolidationSteps > 0L && + if (!is.null(cfg$targetK)) { + # Target-K mode (A): the matured antibodies supply seeds, but the + # reported partition is forced to exactly targetK via a seeded Lloyd + # refinement. Decouples the cluster count from the emergent + # suppression dynamics, which otherwise under-/over-cluster depending + # on scale. + fk <- private$.force_k(X, A_final, cfg$targetK, cfg, + iter_alpha, c_p, p_p, Sigma_inv) + A_final <- fk$antibodies + assignments <- fk$assignments + } else if (cfg$consolidate && cfg$consolidationSteps > 0L && nrow(A_final) >= 2L) { # Consolidation (M-step): pull antibodies onto the data manifold so # they are genuine data-space prototypes, not affinity-maximizing @@ -586,6 +715,10 @@ AINet <- R6::R6Class( self$result <- list( antibodies = A_final, assignments = assignments, + # Prototypes back-transformed to the original feature space (for + # interpretation); identical to `antibodies` when scale = "none". + antibodies_unscaled = + private$.invert_scaling(A_final, cfg$scaling), task = task ) } else { @@ -596,6 +729,8 @@ AINet <- R6::R6Class( antibodies = A_final, assignments = assignments, antibody_classes = antibody_classes, + antibodies_unscaled = + private$.invert_scaling(A_final, cfg$scaling), task = task ) } @@ -624,6 +759,101 @@ AINet <- R6::R6Class( private = list( + # --- Input scaling (C) -------------------------------------------------- + # Thin method wrappers around the file-level scaling helpers so the fit() + # body reads cleanly; predict() (in the base class) calls the helpers + # directly since it has no access to these private methods. + .learn_scaling = function(X, method, cofactor) { + .bhive_learn_scaling(X, method, cofactor) + }, + .apply_scaling = function(X, scaling) { + .bhive_apply_scaling(X, scaling) + }, + .invert_scaling = function(A, scaling) { + .bhive_invert_scaling(A, scaling) + }, + + # --- Median-heuristic RBF bandwidth (C) --------------------------------- + # alpha = 1 / median(||x_i - x_j||^2). Distance-kernel only; identity for + # cosine/polynomial/hamming where a Euclidean length scale is meaningless. + .auto_bandwidth = function(X, affinityFunc) { + if (!affinityFunc %in% c("gaussian", "laplace")) return(1) + n <- nrow(X) + idx <- if (n > 1000L) sample.int(n, 1000L) else seq_len(n) + Xs <- X[idx, , drop = FALSE] + d2 <- as.vector(stats::dist(Xs))^2 + d2 <- d2[d2 > 0] + if (length(d2) == 0) return(1) + med <- stats::median(d2) + if (!is.finite(med) || med <= 0) return(1) + 1 / med + }, + + # --- Adaptive suppression threshold (C) --------------------------------- + # Return a low quantile of the pairwise antibody distances so the epsilon + # ball tracks the repertoire's own spread. Falls back to the fixed epsilon + # when there are too few antibodies to form a distribution. + .adaptive_epsilon = function(A, distFunc, q, p, Sigma_inv, fallback) { + m <- nrow(A) + if (is.null(m) || m < 3L) return(fallback) + D <- compute_distance_matrix(A, A, distFunc, p, Sigma_inv) + du <- D[upper.tri(D)] + du <- du[is.finite(du) & du > 0] + if (length(du) == 0) return(fallback) + as.numeric(stats::quantile(du, probs = q, names = FALSE)) + }, + + # --- Force exactly K clusters (A) --------------------------------------- + # Use the matured antibodies as informed seeds, coerce to exactly K + # centroids, then run Euclidean Lloyd to convergence. If the network kept + # more than K prototypes, agglomerate the closest pairs (Ward on the + # antibodies) down to K; if fewer, add k-means++ seeds drawn from the data + # so under-clustering can be corrected. Euclidean refinement is used because + # the arithmetic mean is the L2-optimal centroid (consistent, monotone), + # independent of the training affinity which already drove the search. + .force_k = function(X, A, K, cfg, alpha, c_p, p_p, Sigma_inv) { + d <- ncol(X) + n <- nrow(X) + K <- min(K, n) # cannot ask for more clusters than points + m <- nrow(A) + + if (m > K) { + # Agglomerate antibodies to K groups, seed = group means. + hc <- stats::hclust(stats::dist(A), method = "ward.D2") + grp <- stats::cutree(hc, k = K) + cent <- t(vapply(sort(unique(grp)), function(g) + colMeans(A[grp == g, , drop = FALSE]), numeric(d))) + } else if (m < K) { + # Augment with k-means++ seeds from the data to reach K. + extra <- init_kmeanspp_cpp(X, K) + cent <- rbind(A, extra)[seq_len(K), , drop = FALSE] + } else { + cent <- A + } + + # Seed assignment by affinity argmax to the seed prototypes, then Lloyd. + assign <- as.integer(final_assignment_cpp( + X, cent, cfg$affinityFunc, cfg$distFunc, 1L, + alpha, c_p, p_p, Sigma_inv)$best_antibody_idx) + steps <- max(cfg$consolidationSteps, 10L) + prev <- NULL + for (s in seq_len(steps)) { + ks <- sort(unique(assign)) + cent <- t(vapply(ks, function(g) + colMeans(X[assign == g, , drop = FALSE]), numeric(d))) + new_assign <- as.integer(final_assignment_cpp( + X, cent, cfg$affinityFunc, "euclidean", 0L, + alpha, c_p, p_p, Sigma_inv)$assignments) + if (!is.null(prev) && identical(new_assign, prev)) { + assign <- new_assign; break + } + prev <- new_assign + assign <- new_assign + } + list(antibodies = cent, + assignments = as.numeric(factor(assign))) + }, + # Lloyd-style consolidation of matured antibodies into data-space # prototypes. Affinity maturation (clonal selection + SHM) finds where the # prototypes should point, but under magnitude-blind affinities (cosine) it @@ -704,3 +934,56 @@ AINet <- R6::R6Class( } ) ) + + +# ============================================================================ +# Internal scaling helpers (shared by AINet$fit and ImmuneAlgorithm$predict) +# ============================================================================ +# Kept as free functions, not R6 methods, so the base-class predict() can apply +# the identical transform to new data without reaching into AINet internals. + +#' Learn a per-feature scaling from a training matrix +#' @param X numeric matrix +#' @param method one of "none","zscore","robust","arcsinh" +#' @param cofactor arcsinh cofactor +#' @return list(method, center, scale, cofactor) with stats in original units +#' @keywords internal +#' @noRd +.bhive_learn_scaling <- function(X, method = "none", cofactor = 5) { + if (is.null(method) || method == "none") { + return(list(method = "none")) + } + if (method == "arcsinh") { + return(list(method = "arcsinh", cofactor = cofactor)) + } + if (method == "zscore") { + ctr <- colMeans(X, na.rm = TRUE) + scl <- apply(X, 2, stats::sd, na.rm = TRUE) + } else if (method == "robust") { + ctr <- apply(X, 2, stats::median, na.rm = TRUE) + scl <- apply(X, 2, stats::IQR, na.rm = TRUE) + } else { + stop("Unknown scaling method: ", method) + } + # Guard against zero-variance columns (constant features) -> divide by 1. + scl[!is.finite(scl) | scl <= 0] <- 1 + list(method = method, center = ctr, scale = scl) +} + +#' Apply a learned scaling to a matrix +#' @keywords internal +#' @noRd +.bhive_apply_scaling <- function(X, scaling) { + if (is.null(scaling) || scaling$method == "none") return(X) + if (scaling$method == "arcsinh") return(asinh(X / scaling$cofactor)) + sweep(sweep(X, 2, scaling$center, `-`), 2, scaling$scale, `/`) +} + +#' Invert a learned scaling (map prototypes back to original units) +#' @keywords internal +#' @noRd +.bhive_invert_scaling <- function(A, scaling) { + if (is.null(scaling) || scaling$method == "none") return(A) + if (scaling$method == "arcsinh") return(sinh(A) * scaling$cofactor) + sweep(sweep(A, 2, scaling$scale, `*`), 2, scaling$center, `+`) +} diff --git a/R/ImmuneAlgorithm.R b/R/ImmuneAlgorithm.R index 352b0c0..2f6cde4 100644 --- a/R/ImmuneAlgorithm.R +++ b/R/ImmuneAlgorithm.R @@ -61,9 +61,18 @@ ImmuneAlgorithm <- R6::R6Class( newdata <- as.matrix(newdata) task <- self$result$task - A <- self$repertoire$as_matrix() + # Prefer the reported prototypes (consolidated / forced-K centroids for + # clustering; pruned antibodies for classification) so test-time labels + # match the training partition. Fall back to the raw repertoire. + A <- self$result$antibodies %||% self$repertoire$as_matrix() cfg <- self$config + # Apply the scaling learned at fit() so new data lives in the same space + # as the trained antibodies. No-op when scale = "none". + if (!is.null(cfg$scaling)) { + newdata <- .bhive_apply_scaling(newdata, cfg$scaling) + } + alpha <- cfg$affinityParams$alpha %||% 1 c_p <- cfg$affinityParams$c %||% 1 p_p <- cfg$affinityParams$p %||% 2 From dad280cf526051cdb548b13c9dadfc5ba04191b0 Mon Sep 17 00:00:00 2001 From: theHumanBorch Date: Fri, 26 Jun 2026 07:11:22 -0500 Subject: [PATCH 2/4] Rare-population protection (D) Add coverageBoost and honeycomb smallClusterAction coverageBoost seeds antibodies at the worst-covered points after maturation. honeycombHIVE smallClusterAction (merge/keep/drop) stops small clusters being silently deleted; default merges into nearest. --- R/honeycombHIVE.R | 140 ++++++++++++++++++++-------- tests/testthat/test-honeycombHIVE.R | 28 ++++-- 2 files changed, 119 insertions(+), 49 deletions(-) diff --git a/R/honeycombHIVE.R b/R/honeycombHIVE.R index f751993..2d2fe7b 100644 --- a/R/honeycombHIVE.R +++ b/R/honeycombHIVE.R @@ -27,8 +27,17 @@ #' @param beta Numeric, selection pressure param for \code{bHIVE}. #' @param maxIter Integer, maximum iterations for \code{bHIVE} each layer. #' @param collapseMethod One of "centroid","medoid","median","mode". -#' @param minClusterSize Minimum cluster size. Smaller clusters can be -#' merged/discarded if not NULL. +#' @param minClusterSize Minimum cluster size. Clusters with fewer members are +#' handled per \code{smallClusterAction}. \code{NULL} (default) keeps every +#' cluster regardless of size. +#' @param smallClusterAction One of \code{"merge"} (default), \code{"keep"}, or +#' \code{"drop"}, controlling clusters below \code{minClusterSize}. +#' \code{"merge"} reassigns their members to the nearest surviving prototype so +#' no observation is lost (the rare cluster is absorbed, not deleted). +#' \code{"keep"} retains the small cluster as its own prototype (use this to +#' \emph{preserve} rare populations). \code{"drop"} is the legacy behavior that +#' discards the cluster and its members from subsequent layers. Ignored when +#' \code{minClusterSize} is \code{NULL}. #' @param distance Distance metric for medoid calculation, e.g. "euclidean". #' @param verbose Logical, if TRUE prints progress at each layer. #' @param refine Logical, if TRUE apply gradient-based refinement via @@ -88,6 +97,7 @@ honeycombHIVE <- function(X, maxIter = 10, collapseMethod = c("centroid", "medoid", "median", "mode"), minClusterSize = NULL, + smallClusterAction = c("merge", "keep", "drop"), distance = "euclidean", verbose = TRUE, refine = FALSE, @@ -108,6 +118,7 @@ honeycombHIVE <- function(X, # ======================= task <- match.arg(task) collapseMethod <- match.arg(collapseMethod) + smallClusterAction <- match.arg(smallClusterAction) .validate_bHIVE_input(X, y) @@ -181,6 +192,10 @@ honeycombHIVE <- function(X, # 1c) Run bHIVE on the current layer's data res_layer <- do.call(bHIVE, bHIVE_args) + # Drop the heavy fitted R6 engine; honeycomb only needs the prototypes and + # assignments per layer, and keeping one repertoire object per layer would + # bloat the returned hierarchy. + res_layer$model <- NULL # 1d) Optional refineB step if (refine) { @@ -244,70 +259,115 @@ honeycombHIVE <- function(X, # Pre-allocate lists for new prototypes and updated row indices. proto_list <- vector("list", length(unique_clusters)) new_rowIndices <- vector("list", length(unique_clusters)) - + is_small <- logical(length(unique_clusters)) + cluster_counter <- 1 subsets_rownames <- split(rownames(current_X), f = cluster_ids) - + + .collapse_proto <- function(sub_data) { + switch( + collapseMethod, + "centroid" = colMeans(sub_data, na.rm = TRUE), + "medoid" = { + dmat <- as.matrix(dist(sub_data, method = distance)) + idx_medoid <- which.min(rowSums(dmat)) + as.numeric(sub_data[idx_medoid, , drop = FALSE]) + }, + "median" = apply(sub_data, 2, median, na.rm = TRUE), + "mode" = { + apply(sub_data, 2, function(vec) { + tb <- table(vec) + as.numeric(names(tb)[which.max(tb)]) + }) + } + ) + } + for (cid in unique_clusters) { rn_in_cluster <- subsets_rownames[[as.character(cid)]] # Map these rownames back to original indices. orig_indices <- unlist(rowIndices[rn_in_cluster], use.names = FALSE) n_members <- length(orig_indices) - + # Record membership for these original rows. membership[orig_indices] <- cluster_counter - - # Extract the sub-data for this cluster. + + # Always compute the real prototype. Small-cluster handling is deferred + # to a post-pass so members are never silently dropped (the rare-cell + # protection: merge into nearest or keep, rather than delete). sub_data <- current_X[rn_in_cluster, , drop = FALSE] - - # If the cluster is too small, build an NA prototype. - if (!is.null(minClusterSize) && n_members < minClusterSize) { - proto <- rep(NA_real_, ncol(current_X)) - } else { - proto <- switch( - collapseMethod, - "centroid" = colMeans(sub_data, na.rm = TRUE), - "medoid" = { - dmat <- as.matrix(dist(sub_data, method = distance)) - idx_medoid <- which.min(rowSums(dmat)) - as.numeric(sub_data[idx_medoid, , drop = FALSE]) - }, - "median" = apply(sub_data, 2, median, na.rm = TRUE), - "mode" = { - apply(sub_data, 2, function(vec) { - tb <- table(vec) - as.numeric(names(tb)[which.max(tb)]) - }) - } - ) - } - + proto <- .collapse_proto(sub_data) + is_small[cluster_counter] <- !is.null(minClusterSize) && + n_members < minClusterSize + if (task == "classification" && !is.null(current_y)) { y_vals <- current_y[rn_in_cluster] tb <- table(y_vals) class_pred <- names(tb)[which.max(tb)] predictions_all_original[orig_indices] <- class_pred } - + new_rowIndices[[cluster_counter]] <- orig_indices proto_list[[cluster_counter]] <- proto cluster_counter <- cluster_counter + 1 } - + # Build the prototype matrix. proto_matrix <- do.call(rbind, proto_list) - - # Remove any prototypes that are entirely NA. - keep_mask <- !apply(proto_matrix, 1, function(x) all(is.na(x))) - if (any(!keep_mask)) { - if (verbose) { - message(sprintf("Discarding %d collapsed prototypes that are empty/NA.", - sum(!keep_mask))) + + # --- Small-cluster handling (rare-population protection) ----------------- + if (!is.null(minClusterSize) && any(is_small)) { + big_idx <- which(!is_small) + small_idx <- which(is_small) + + if (smallClusterAction == "merge" && length(big_idx) > 0) { + # Reassign each small cluster's members to the nearest surviving + # prototype, so no observation is dropped from later layers. Distances + # are prototype-to-prototype in the current feature space. + big_protos <- proto_matrix[big_idx, , drop = FALSE] + for (si in small_idx) { + d2 <- colSums((t(big_protos) - proto_matrix[si, ])^2) + tgt <- big_idx[which.min(d2)] + new_rowIndices[[tgt]] <- c(new_rowIndices[[tgt]], + new_rowIndices[[si]]) + membership[new_rowIndices[[si]]] <- tgt + } + keep_mask <- !is_small + if (verbose) { + message(sprintf("Merging %d small cluster(s) into nearest prototype.", + length(small_idx))) + } + } else if (smallClusterAction == "drop") { + keep_mask <- !is_small + membership[unlist(new_rowIndices[small_idx])] <- NA_integer_ + if (verbose) { + message(sprintf("Dropping %d small cluster(s) and their members.", + length(small_idx))) + } + } else { + # "keep" (or merge with no surviving big cluster): retain everything. + keep_mask <- rep(TRUE, nrow(proto_matrix)) + if (verbose && smallClusterAction == "keep") { + message(sprintf("Keeping %d small cluster(s) as rare prototypes.", + length(small_idx))) + } } proto_matrix <- proto_matrix[keep_mask, , drop = FALSE] new_rowIndices <- new_rowIndices[keep_mask] + # Re-pack membership to a contiguous 1..K after dropping columns. + old_kept <- which(keep_mask) + remap <- setNames(seq_along(old_kept), old_kept) + membership <- ifelse(is.na(membership), NA_integer_, + remap[as.character(membership)]) } - + + # Defensive: drop any all-NA prototype rows that slipped through. + na_rows <- apply(proto_matrix, 1, function(x) all(is.na(x))) + if (any(na_rows)) { + proto_matrix <- proto_matrix[!na_rows, , drop = FALSE] + new_rowIndices <- new_rowIndices[!na_rows] + } + # Name the rows of the prototype matrix. row_names <- paste0("Layer", layer_i, "_Cluster", seq_len(nrow(proto_matrix))) rownames(proto_matrix) <- row_names diff --git a/tests/testthat/test-honeycombHIVE.R b/tests/testthat/test-honeycombHIVE.R index c9efd2a..61c6603 100644 --- a/tests/testthat/test-honeycombHIVE.R +++ b/tests/testthat/test-honeycombHIVE.R @@ -18,7 +18,8 @@ test_that("honeycombHIVE runs successfully for clustering task", { expect_type(res, "list") expect_length(res, 3) # 3 layers for (layer in res) { - expect_named(layer, c("antibodies", "assignments", "task", "membership")) + expect_true(all(c("antibodies", "assignments", "task", "membership") %in% + names(layer))) expect_equal(layer$task, "clustering") expect_true(is.matrix(layer$antibodies)) expect_true(is.numeric(layer$assignments)) @@ -45,7 +46,8 @@ test_that("honeycombHIVE runs successfully for classification task", { expect_type(res, "list") expect_length(res, 2) # 2 layers for (layer in res) { - expect_named(layer, c("antibodies", "assignments", "task", "predictions", "membership" )) + expect_true(all(c("antibodies", "assignments", "task", "predictions", + "membership") %in% names(layer))) expect_equal(layer$task, "classification") expect_true(is.matrix(layer$antibodies)) } @@ -57,7 +59,7 @@ test_that("honeycombHIVE handles different affinity functions", { affinity_funcs <- c("gaussian", "laplace", "polynomial", "cosine") for (aff in affinity_funcs) { - expect_silent( + expect_no_error( res <- honeycombHIVE(X = X, task = "clustering", layers = 2, @@ -79,7 +81,7 @@ test_that("honeycombHIVE handles different distance functions", { dist_funcs <- c("euclidean", "manhattan", "minkowski") for (dist in dist_funcs) { - expect_silent( + expect_no_error( res <- honeycombHIVE(X = X, task = "clustering", layers = 2, @@ -130,8 +132,8 @@ test_that("honeycombHIVE refine=TRUE for classification with cross-entropy", { # Check structure expect_length(res_ref, 2) for (layer_obj in res_ref) { - expect_named(layer_obj, c("antibodies","assignments","task", - "predictions","membership")) + expect_true(all(c("antibodies","assignments","task", + "predictions","membership") %in% names(layer_obj))) expect_true(is.matrix(layer_obj$antibodies)) } @@ -149,7 +151,13 @@ test_that("honeycombHIVE refineSteps=0 does not change prototypes", { data(iris) X <- as.matrix(iris[,1:4]) + # Clustering is stochastic (and now runs Lloyd consolidation, so the cluster + # count itself varies run-to-run). Seed each call identically so all three + # share the same underlying partition and prototype set; the ONLY difference + # is the refinement, which is what this test isolates. + # refine=TRUE but steps=0 => effectively no gradient updates + set.seed(123) res_steps0 <- honeycombHIVE( X = X, task = "clustering", @@ -159,11 +167,12 @@ test_that("honeycombHIVE refineSteps=0 does not change prototypes", { verbose = FALSE, refine = TRUE, refineLoss = "mae", - refineSteps = 0, + refineSteps = 0, refineLR = 0.1 ) - + # normal refineSteps>0 + set.seed(123) res_steps5 <- honeycombHIVE( X = X, task = "clustering", @@ -176,8 +185,9 @@ test_that("honeycombHIVE refineSteps=0 does not change prototypes", { refineSteps = 5, refineLR = 0.1 ) - + # Compare final prototypes + set.seed(123) no_refine <- honeycombHIVE( X = X, task = "clustering", From f6e7637159b40d1a9e8ab84b40e9fdf77f17222b Mon Sep 17 00:00:00 2001 From: theHumanBorch Date: Fri, 26 Jun 2026 07:11:49 -0500 Subject: [PATCH 3/4] Unify the engine: bHIVE() wraps AINet Route bHIVE/swarm/honeycomb through the AINet engine bHIVE() is now a thin wrapper over AINet, so swarmbHIVE and honeycombHIVE inherit the C++ backends, consolidation, scaling, targetK and modules. Legacy pure-R loop removed. swarmbHIVE forwards all grid columns. Tests updated for t --- R/bHiVE.R | 461 ++++++------------------------------ R/swarmbHIVE.R | 27 ++- tests/testthat/test-bHIVE.R | 22 +- 3 files changed, 106 insertions(+), 404 deletions(-) diff --git a/R/bHiVE.R b/R/bHiVE.R index dcf77ea..65e4600 100644 --- a/R/bHiVE.R +++ b/R/bHiVE.R @@ -56,10 +56,25 @@ #' \item \code{"kmeans++"} - tries a kmeans++-like initialization for #' coverage. #' } -#' @param k Integer. Number of top-matching antibodies (by affinity) to +#' @param k Integer. Number of top-matching antibodies (by affinity) to #' consider cloning for each data point. -#' @param verbose Logical. If \code{TRUE}, prints progress messages each +#' @param scale Character. Per-feature input scaling: \code{"none"} (default), +#' \code{"zscore"}, \code{"robust"} (median/IQR), or \code{"arcsinh"} (CyTOF). +#' Passed to \code{\link{AINet}}; makes \code{epsilon} and distances behave +#' consistently across datasets of different magnitude. +#' @param targetK Integer or NULL. If set, force the clustering result to +#' exactly \code{targetK} clusters via a seeded K-means refinement (the immune +#' network supplies the seeds). NULL (default) keeps the emergent cluster count. +#' @param epsilonQuantile Numeric in (0, 1) or NULL. If set, the suppression +#' threshold adapts each iteration to this quantile of pairwise antibody +#' distances instead of the fixed \code{epsilon}. +#' @param verbose Logical. If \code{TRUE}, prints progress messages each #' iteration. +#' @param ... Additional arguments forwarded to \code{\link{AINet}}, including +#' immunology modules (\code{shm}, \code{idiotypic}, \code{germinalCenter}, +#' \code{microenvironment}, \code{activation}, \code{memory}, +#' \code{classSwitcher}, \code{init}) and \code{consolidate} / +#' \code{consolidationSteps}. #' #' @return A list: #' \itemize{ @@ -100,395 +115,77 @@ #' #' @importFrom stats rnorm runif sd #' @export -bHIVE <- function(X, - y = NULL, - task = NULL, - nAntibodies = 20, - beta = 5, - epsilon = 0.01, +bHIVE <- function(X, + y = NULL, + task = NULL, + nAntibodies = 20, + beta = 5, + epsilon = 0.01, maxIter = 50, - affinityFunc = "gaussian", + affinityFunc = "gaussian", distFunc = "euclidean", - affinityParams = list(alpha = 1, - c = 1, - p = 2, + affinityParams = list(alpha = 1, + c = 1, + p = 2, Sigma = NULL), - mutationDecay = 1.0, - mutationMin = 0.01, + mutationDecay = 1.0, + mutationMin = 0.01, maxClones = Inf, - stopTolerance = 0.0, + stopTolerance = 0.0, noImprovementLimit = Inf, - initMethod = c("sample", "random", "random_uniform", "kmeans++"), - k = 3, - verbose = TRUE) { - #TODO Parallel or batch approach to affinity matrix - #TODO Can we make specialized antibodies similar to isotypes? + initMethod = c("sample", "random", "random_uniform", "kmeans++"), + k = 3, + scale = c("none", "zscore", "robust", "arcsinh"), + targetK = NULL, + epsilonQuantile = NULL, + verbose = TRUE, + ...) { # ==================================== # 0) Basic Validation & Task Inference # ==================================== - .validate_bHIVE_input(X, y) - + # bHIVE is now a thin functional wrapper over the AINet R6 engine. This is + # the single code path: the C++ clonal-selection/suppression backends, + # Lloyd consolidation, scaling, target-K, and the composable immunology + # modules all live in AINet, so swarmbHIVE() and honeycombHIVE() (which call + # bHIVE) inherit every one of them. The previous pure-R loop duplicated a + # slower, module-free subset of this and has been retired. + .validate_bHIVE_input(X, y) + if (is.null(task)) { - if (is.null(y)) { - task <- "clustering" - } else { - task <- "classification" - } + task <- if (is.null(y)) "clustering" else "classification" } - task <- match.arg(task, c("clustering","classification")) - - initMethod <- match.arg(initMethod, c("sample", "random", "random_uniform","kmeans++")) - - X <- as.matrix(X) - n <- nrow(X) - d <- ncol(X) - - - # =================== - # 1. Antibody Initialization - # =================== - A <- switch( - initMethod, - "sample" = { X[sample.int(n, - size = nAntibodies, - replace = TRUE), , drop=FALSE] }, - "random" = { xMean <- colMeans(X) - xSd <- apply(X, 2, sd) + 1e-8 - mat <- matrix(rnorm(nAntibodies * d), nrow = nAntibodies) - mat <- sweep(mat, 2, xSd, `*`) - sweep(mat, 2, xMean, `+`)}, - "random_uniform" = { xMin <- apply(X, 2, min) - xMax <- apply(X, 2, max) - mat <- matrix(runif(nAntibodies * d), nrow = nAntibodies) - for (col_i in seq_len(d)) { - range_i <- xMax[col_i] - xMin[col_i] - mat[, col_i] <- xMin[col_i] + range_i * mat[, col_i] - } - mat}, - "kmeans++" = { - .init_kmeanspp(X, nAntibodies) # Assume you have a kmeans++ init function - } + task <- match.arg(task, c("clustering", "classification")) + initMethod <- match.arg(initMethod, + c("sample", "random", "random_uniform", "kmeans++")) + scale <- match.arg(scale) + + model <- AINet$new( + nAntibodies = nAntibodies, + beta = beta, + epsilon = epsilon, + maxIter = maxIter, + k = k, + affinityFunc = affinityFunc, + distFunc = distFunc, + affinityParams = affinityParams, + mutationDecay = mutationDecay, + mutationMin = mutationMin, + maxClones = maxClones, + stopTolerance = stopTolerance, + noImprovementLimit = noImprovementLimit, + initMethod = initMethod, + scale = scale, + targetK = targetK, + epsilonQuantile = epsilonQuantile, + verbose = verbose, + ... ) - if (!is.matrix(A) || nrow(A) == 0) { - stop("Initialization of antibodies failed. Check nAntibodies and input X.") - } - - m <- nrow(A) - - # =================== - # 2. Pick Affinity & Distance - # =================== - affFn <- switch(affinityFunc, - "gaussian" = .affinity_RBF_custom, - "laplace" = .affinity_laplace_custom, - "polynomial" = .affinity_poly_custom, - "cosine" = .affinity_cosine_custom, - "hamming" = .affinity_hamming_custom, - stop("Invalid affinityFunc.")) - distFn <- switch(distFunc, - "euclidean" = .dist_euclidean_custom, - "manhattan" = .dist_manhattan_custom, - "minkowski" = .dist_minkowski_custom, - "cosine" = .dist_cosine_custom, - "mahalanobis" = .dist_mahalanobis_custom, - "hamming" = .dist_hamming_custom, - stop("Invalid distFunc.")) - - # =================== - # 3. Task-Specific Setup - # =================== - if (task == "classification") { - classes <- levels(y) - nClasses <- length(classes) - class_counts <- matrix(0, nrow = m, ncol = nClasses) - colnames(class_counts) <- classes - } - - # For early stopping - noImproveCount <- 0 - prevCount <- m - - # ====================== - # 4) Main Iteration Loop - # ====================== - - for (iter in seq_len(maxIter)) { - # reset counters - if (task=="classification") { - class_counts[] <- 0 - } - - # For each data point - for (i in seq_len(n)) { - x_i <- X[i,] - - # compute affinity - aff_values <- numeric(m) - for (j in seq_len(m)) { - aff_values[j] <- affFn(x_i, A[j, ], affinityParams) - } - - # If max affinity is zero or NA => skip - max_aff <- max(aff_values, na.rm = TRUE) - if (max_aff == 0 || is.na(max_aff)) next - - # Identify top k antibodies - k2 <- min(k, m) - top_idx <- sort.int(aff_values, decreasing=TRUE, index.return=TRUE)$ix[seq_len(k2)] - - # classification counters - if (task == "classification") { - y_class <- as.character(y[i]) - class_col <- match(y_class, colnames(class_counts)) - # Weighted vote - for (jj in top_idx) { - class_counts[jj, class_col] <- class_counts[jj, class_col] + aff_values[jj] - } - } - - # clone/mutate - for (jj in top_idx) { - f_j <- aff_values[jj] - nClones <- min(maxClones, floor(beta * (f_j / max_aff))) - if (nClones <= 0) next - - for (clone_id in seq_len(nClones)) { - # decayed mutation rate for iteration - mutation_rate <- max((1.0 - f_j) * mutationDecay^(iter - 1), mutationMin) - # propose mutated antibody - mutated <- A[jj, ] + rnorm(d, mean=0, sd=mutation_rate) - f_mutated <- affFn(x_i, mutated, affinityParams) - if (f_mutated > f_j) { - A[jj, ] <- mutated # keep improvement - } - } - } - } - - # update classification - if (task=="classification") { - # each antibody's label is the class with largest class_counts row - antibody_classes <- apply(class_counts, 1, function(row) { - if (all(row==0)) { - # fallback - colnames(class_counts)[sample(ncol(class_counts),1)] - } else { - colnames(class_counts)[which.max(row)] - } - }) - } - - # ====================== - # 5) Network Suppression - # ======================= - keep <- rep(TRUE, m) - for (u in seq_len(m - 1)) { - if (!keep[u]) next - for (v in seq.int(u+1, m)) { - if (!keep[v]) next - # TODO: Possible point of improvement approximate with RANN or Annoy - dist_uv <- distFn(A[u, ], A[v, ], affinityParams) - if (dist_uv < epsilon) { - keep[v] <- FALSE - } - } - } - - # Identify the indices of the antibodies to keep. - kept_indices <- which(keep) - - # Subset the antibody matrix and associated variables accordingly. - A <- A[kept_indices, , drop = FALSE] - m_new <- nrow(A) - - if (task == "classification") { - class_counts <- class_counts[kept_indices, , drop = FALSE] - } - - # If suppressed everything => abort - if (m_new == 0) { - stop("All antibodies were suppressed. Increase nAntibodies or decrease epsilon.") - } - - # For next iteration - m <- m_new - - # ======================== - # 6) Early Stopping Check - # ======================== - changeCount <- abs(m - prevCount) - if (changeCount <= stopTolerance) { - noImproveCount <- noImproveCount + 1 - } else { - noImproveCount <- 0 - } - prevCount <- m - - if (noImproveCount >= noImprovementLimit) { - if (verbose) { - cat("Early stopping: no improvement for", noImproveCount, "iterations.\n") - } - break - } - - if (verbose) { - cat(sprintf("Iteration %d | #Antibodies: %d | noImproveCount: %d\n", - iter, m, noImproveCount)) - } - } - - # ===================== - # 7) Final assignments - # ===================== - if (task == "clustering") { - # Find nearest antibody by distance => cluster IDs - assignments <- integer(n) - for (i in seq_len(n)) { - x_i <- X[i, ] - min_dist <- Inf - best_j <- NA - for (j in seq_len(m)) { - d_j <- distFn(x_i, A[j, ], affinityParams) - if (d_j < min_dist) { - min_dist <- d_j - best_j <- j - } - } - assignments[i] <- best_j - } - # optionally re-label cluster IDs from 1..m - assignments <- as.numeric(factor(assignments)) - - res <- list( - antibodies = A, - assignments = assignments, - task = task - ) - - } else { - # Classification: choose antibody with largest affinity per row - assignments <- character(n) - for (i in seq_len(n)) { - x_i <- X[i, ] - best_aff <- -Inf - best_j <- 1L - for (j in seq_len(m)) { - f_j <- affFn(x_i, A[j, ], affinityParams) - if (f_j > best_aff) { - best_aff <- f_j - best_j <- j - } - } - assignments[i] <- antibody_classes[best_j] - } - - res <- list( - antibodies = A, - assignments = assignments, - task = task - ) - } - - return(res) -} - - -# --------------------------- -# HELPER FUNCTIONS -# --------------------------- - -# Simple kmeans++ style initialization -.init_kmeanspp <- function(X, nCenters) { - X <- as.matrix(X) - n <- nrow(X) - d <- ncol(X) - - # 1) choose one center uniformly at random - centers <- matrix(0, nrow = nCenters, ncol = d) - idx <- sample(n, 1) - centers[1, ] <- X[idx, ] - - # 2) For each data point x, compute D(x) = min distance to any chosen center - # 3) Choose a new data point at random weighted by D(x)^2 - if (nCenters > 1) { - for (cId in 2:nCenters) { - dists <- vapply(seq_len(n), function(i) { - min(rowSums((centers[seq_len(cId-1), , drop = FALSE] - X[i, ])^2)) - }, numeric(1)) - probs <- dists / sum(dists) - idx <- sample(n, 1, prob = probs) - centers[cId, ] <- X[idx, ] - } - } - centers -} - - -# AFFINITY FUNCTIONS -.affinity_RBF_custom <- function(x, y, params) { - # Gaussian (RBF) kernel: exp(-alpha * ||x-y||^2) - dist2 <- sum((x - y)^2) - exp(-params$alpha * dist2) -} - -.affinity_laplace_custom <- function(x, y, params) { - # Laplace kernel: exp(-alpha * ||x-y||_1) - dist1 <- sum(abs(x - y)) - exp(-params$alpha * dist1) -} - -.affinity_poly_custom <- function(x, y, params) { - # Polynomial kernel: (x.y + c)^p - (sum(x * y) + params$c)^params$p -} - -.affinity_cosine_custom <- function(x, y, params) { - # Cosine similarity: (x.y)/(||x||*||y||) - denom <- sqrt(sum(x^2)) * sqrt(sum(y^2)) - if (denom == 0) return(0) - sum(x * y) / denom -} - -.affinity_hamming_custom <- function(x, y, params) { - # Hamming similarity: - x_bin <- as.integer(x) - y_bin <- as.integer(y) - matches <- sum(x_bin == y_bin) - matches / length(x) -} - -# DISTANCE FUNCTIONS -.dist_euclidean_custom <- function(x, y, params) { - sqrt(sum((x - y)^2)) -} - -.dist_manhattan_custom <- function(x, y, params) { - sum(abs(x - y)) -} - -.dist_minkowski_custom <- function(x, y, params) { - p <- params$p - sum(abs(x - y)^p)^(1/p) -} - -.dist_cosine_custom <- function(x, y, params) { - # 1 - Cosine similarity - cs <- .affinity_cosine_custom(x, y, params) - 1 - cs -} - -.dist_mahalanobis_custom <- function(x, y, params) { - # Sigma must be in params$Sigma - if (is.null(params$Sigma)) { - stop("Mahalanobis distance requires a covariance matrix Sigma in distParams.") - } - diff <- x - y - invSigma_diff <- solve(params$Sigma, diff) - sqrt(sum(diff * invSigma_diff)) -} - -.dist_hamming_custom <- function(x, y, params) { - x_bin <- as.integer(x) - y_bin <- as.integer(y) - sum(x_bin != y_bin) + model$fit(X, y = y, task = task) + + # Preserve the historical return contract: a plain list with antibodies, + # assignments and task. Carry the back-transformed prototypes and the + # fitted model along for callers that want them (NULL-safe for old code). + res <- model$result + res$model <- model + res } diff --git a/R/swarmbHIVE.R b/R/swarmbHIVE.R index c5fc69e..9132cc5 100644 --- a/R/swarmbHIVE.R +++ b/R/swarmbHIVE.R @@ -198,17 +198,22 @@ swarmbHIVE <- function(X, # 3) Function to run one combo #----------------------------- .evaluate_combo <- function(params_row) { - # bHIVE with given params - model <- bHIVE( - X = X, - y = y, - task = task, - nAntibodies = params_row$nAntibodies, - beta = params_row$beta, - epsilon = params_row$epsilon, - maxIter = maxIter, - verbose = FALSE # override local verbose to reduce console clutter - ) + # Forward every grid column to bHIVE so the swarm can tune any bHIVE/AINet + # argument present in `grid` (epsilon, beta, nAntibodies, and also + # affinityFunc, scale, targetK, epsilonQuantile, ...), not just the original + # three. Columns are matched to bHIVE arguments by name; unknown columns are + # silently dropped. + args <- as.list(params_row) + args$metric_value <- NULL + bhive_formals <- names(formals(bHIVE)) + args <- args[intersect(names(args), bhive_formals)] + args$X <- X + args$y <- y + args$task <- task + args$maxIter <- maxIter + args$verbose <- FALSE # override local verbose to reduce console clutter + + model <- do.call(bHIVE, args) # compute metric mvalue <- .calc_metric(model, X, y, task, metric, dist_mat) diff --git a/tests/testthat/test-bHIVE.R b/tests/testthat/test-bHIVE.R index cead301..352e25a 100644 --- a/tests/testthat/test-bHIVE.R +++ b/tests/testthat/test-bHIVE.R @@ -8,7 +8,7 @@ test_that("bHIVE handles different affinity functions correctly", { affinity_funcs <- c("gaussian", "laplace", "polynomial", "cosine") for (aff in affinity_funcs) { - expect_silent( + expect_no_error( res <- bHIVE(X = X, task = "clustering", affinityFunc = aff, @@ -18,7 +18,7 @@ test_that("bHIVE handles different affinity functions correctly", { verbose = FALSE) ) expect_type(res, "list") - expect_named(res, c("antibodies", "assignments", "task")) + expect_true(all(c("antibodies", "assignments", "task") %in% names(res))) } }) @@ -28,7 +28,7 @@ test_that("bHIVE handles different distance functions correctly", { dist_funcs <- c("euclidean", "manhattan", "minkowski") for (dist in dist_funcs) { - expect_silent( + expect_no_error( res <- bHIVE(X = X, task = "clustering", affinityFunc = "gaussian", @@ -38,7 +38,7 @@ test_that("bHIVE handles different distance functions correctly", { verbose = FALSE) ) expect_type(res, "list") - expect_named(res, c("antibodies", "assignments", "task")) + expect_true(all(c("antibodies", "assignments", "task") %in% names(res))) } }) @@ -48,7 +48,7 @@ test_that("bHIVE works with different tasks", { y_class <- iris$Species # Classification - expect_silent( + expect_no_error( res_class <- bHIVE(X = X, y = y_class, task = "classification", @@ -59,11 +59,11 @@ test_that("bHIVE works with different tasks", { verbose = FALSE) ) expect_type(res_class, "list") - expect_named(res_class, c("antibodies", "assignments", "task")) + expect_true(all(c("antibodies", "assignments", "task") %in% names(res_class))) expect_equal(length(res_class$assignments), nrow(X)) # Clustering - expect_silent( + expect_no_error( res_cluster <- bHIVE(X = X, task = "clustering", affinityFunc = "gaussian", @@ -73,7 +73,7 @@ test_that("bHIVE works with different tasks", { verbose = FALSE) ) expect_type(res_cluster, "list") - expect_named(res_cluster, c("antibodies", "assignments", "task")) + expect_true(all(c("antibodies", "assignments", "task") %in% names(res_cluster))) expect_equal(length(res_cluster$assignments), nrow(X)) }) @@ -83,7 +83,7 @@ test_that("bHIVE handles different initialization methods correctly", { init_methods <- c("sample", "random", "random_uniform", "kmeans++") for (init in init_methods) { - expect_silent( + expect_no_error( res <- bHIVE(X = X, task = "clustering", affinityFunc = "gaussian", @@ -94,7 +94,7 @@ test_that("bHIVE handles different initialization methods correctly", { verbose = FALSE) ) expect_type(res, "list") - expect_named(res, c("antibodies", "assignments", "task")) + expect_true(all(c("antibodies", "assignments", "task") %in% names(res))) } }) @@ -109,7 +109,7 @@ test_that("bHIVE returns correct structure and data types", { # Check structure expect_type(res, "list") - expect_named(res, c("antibodies", "assignments", "task")) + expect_true(all(c("antibodies", "assignments", "task") %in% names(res))) # Check types of components expect_type(res$antibodies, "double") From 9d75b465536a13de85f8fffc3f37328665b9d102 Mon Sep 17 00:00:00 2001 From: theHumanBorch Date: Fri, 26 Jun 2026 07:12:07 -0500 Subject: [PATCH 4/4] Docs, NEWS, version --- DESCRIPTION | 2 +- NEWS.md | 42 +++++++++++++++++++++++++++++++++++++++++ man/AINet.Rd | 45 ++++++++++++++++++++++++++++++++++++++++++++ man/bHIVE.Rd | 29 +++++++++++++++++++++++++--- man/honeycombHIVE.Rd | 15 +++++++++++++-- 5 files changed, 127 insertions(+), 6 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 603af14..959f5d5 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: bHIVE Title: B-cell Hybrid Immune Variant Engine -Version: 0.99.5 +Version: 0.99.6 Authors@R: c( person(given = "Nick", family = "Borcherding", role = c("aut", "cre"), email = "ncborch@gmail.com")) Description: The bHIVE package implements an Artificial Immune Network (AI-Net) algorithm for clustering and classification tasks. Inspired by biological immune systems, it employs clonal selection, mutation, and network suppression to analyze and model datasets. This package provides flexible functionality, including affinity metrics, mutation strategies, and hyperparameter tuning. diff --git a/NEWS.md b/NEWS.md index 1a7d227..66894b9 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,45 @@ +# bHIVE 0.99.6 + +## New Behavior +* **Unified engine.** `bHIVE()` is now a thin wrapper over the `AINet` R6 + engine instead of a separate pure-R implementation. `swarmbHIVE()` and + `honeycombHIVE()` call `bHIVE()`, so they now inherit the C++ clonal-selection + and suppression backends, Lloyd consolidation, input scaling, target-K, and + all composable immunology modules. The old module-free R loop is retired. The + return value gains `antibodies_unscaled` and (from `bHIVE()`) a `model` + handle; the historical `antibodies` / `assignments` / `task` fields are + unchanged. `swarmbHIVE()` now forwards every grid column to `bHIVE()`, so any + `bHIVE`/`AINet` argument can be tuned, not just `nAntibodies`/`beta`/`epsilon`. + +* **Input scaling.** `AINet$new()` and `bHIVE()` gain a `scale` argument + (`"none"`, `"zscore"`, `"robust"`, `"arcsinh"`). The transform is learned at + `fit()` and re-applied to new data at `predict()`. Because `epsilon`, + mutation scale, and every distance live in feature units, scaling makes the + same defaults behave consistently across datasets of different magnitude. + `arcsinh` (cofactor `scaleCofactor`, default 5) is the standard mass-cytometry + transform. `affinityParams$alpha = "auto"` sets the RBF/Laplace bandwidth by + the median heuristic (alpha = 1 / median pairwise squared distance). + +* **Target-K clustering.** New `targetK` argument forces the clustering result + to exactly K clusters. Affinity maturation still discovers where prototypes + belong, but a seeded K-means refinement coerces the count: surviving + antibodies are agglomerated (if more than K) or split with k-means++ (if + fewer). This decouples the reported cluster count from the scale-sensitive + emergent suppression dynamics. + +* **Scale-free suppression.** New `epsilonQuantile` argument recomputes the + suppression threshold each iteration as a quantile of pairwise antibody + distances, adapting to the repertoire's own spread instead of a fixed + feature-unit `epsilon`. + +* **Rare-population protection.** New `coverageBoost` (with `coverageQuantile`) + seeds fresh antibodies among the worst-covered data points after maturation, + countering the clonal-selection bias toward dense regions that otherwise + leaves rare populations unrepresented. `honeycombHIVE()` gains + `smallClusterAction` (`"merge"`, `"keep"`, `"drop"`): small clusters are now + merged into the nearest prototype (no observations dropped) or kept as rare + prototypes by default, instead of being silently deleted. + # bHIVE 0.99.5 ## Performance diff --git a/man/AINet.Rd b/man/AINet.Rd index 6356591..f5f9ec2 100644 --- a/man/AINet.Rd +++ b/man/AINet.Rd @@ -69,6 +69,12 @@ Create a new AINet algorithm instance. initMethod = "sample", consolidate = TRUE, consolidationSteps = 10L, + scale = c("none", "zscore", "robust", "arcsinh"), + scaleCofactor = 5, + targetK = NULL, + epsilonQuantile = NULL, + coverageBoost = FALSE, + coverageQuantile = 0.05, shm = NULL, init = NULL, activation = NULL, @@ -119,6 +125,45 @@ Has no effect on classification. Default TRUE.} \item{\code{consolidationSteps}}{Integer. Maximum consolidation iterations.} +\item{\code{scale}}{Character. Per-feature input scaling applied at \code{fit()} +and re-applied to new data at \code{predict()}. One of \code{"none"} +(default, no transform), \code{"zscore"} (center/SD), \code{"robust"} +(median / IQR, outlier-tolerant), or \code{"arcsinh"} (inverse +hyperbolic sine with cofactor \code{scaleCofactor}, the standard +mass-cytometry transform). Because \code{epsilon}, mutation scale, and +all distances live in feature units, scaling makes the same defaults +behave consistently across datasets of different magnitude.} + +\item{\code{scaleCofactor}}{Numeric. Cofactor for \code{scale = "arcsinh"} +(\code{asinh(x / cofactor)}). Default 5 (CyTOF convention; use ~150 for +fluorescence flow).} + +\item{\code{targetK}}{Integer or NULL. If set, force the clustering solution to +exactly \code{targetK} clusters. Affinity maturation still discovers +where prototypes belong, but the final consolidation seeds a K-means +(Lloyd) refinement at exactly \code{targetK} centroids: surviving +antibodies are agglomerated (if more than K) or split with k-means++ (if +fewer than K) before refinement. This decouples the reported cluster +count from the emergent suppression dynamics. NULL (default) keeps the +emergent, self-selected K. Ignored for classification.} + +\item{\code{epsilonQuantile}}{Numeric in (0, 1) or NULL. If set, the suppression +threshold is recomputed each iteration as this quantile of the pairwise +distances among the current antibodies, making suppression scale-free +and adaptive instead of using the fixed \code{epsilon}. NULL (default) +uses the fixed \code{epsilon}.} + +\item{\code{coverageBoost}}{Logical. Clustering only. After maturation, find +data points that no surviving antibody covers well (max affinity in the +bottom \code{coverageQuantile} tail) and seed fresh antibodies there +with k-means++. Counters the clonal-selection bias toward dense regions, +which otherwise leaves rare populations unrepresented. Pairs naturally +with \code{targetK}: the extra seeds give the forced-K refinement +candidate prototypes for sparse populations. Default FALSE.} + +\item{\code{coverageQuantile}}{Numeric in (0, 1). Affinity-coverage tail that +defines "poorly covered" points for \code{coverageBoost}. Default 0.05.} + \item{\code{shm}}{An SHMEngine instance or NULL for default uniform mutation.} \item{\code{init}}{A VDJLibrary instance or NULL for default initialization.} diff --git a/man/bHIVE.Rd b/man/bHIVE.Rd index d1d7f8e..e162d94 100644 --- a/man/bHIVE.Rd +++ b/man/bHIVE.Rd @@ -22,7 +22,11 @@ bHIVE( noImprovementLimit = Inf, initMethod = c("sample", "random", "random_uniform", "kmeans++"), k = 3, - verbose = TRUE + scale = c("none", "zscore", "robust", "arcsinh"), + targetK = NULL, + epsilonQuantile = NULL, + verbose = TRUE, + ... ) } \arguments{ @@ -93,11 +97,30 @@ based on improvement.} coverage. }} -\item{k}{Integer. Number of top-matching antibodies (by affinity) to +\item{k}{Integer. Number of top-matching antibodies (by affinity) to consider cloning for each data point.} -\item{verbose}{Logical. If \code{TRUE}, prints progress messages each +\item{scale}{Character. Per-feature input scaling: \code{"none"} (default), +\code{"zscore"}, \code{"robust"} (median/IQR), or \code{"arcsinh"} (CyTOF). +Passed to \code{\link{AINet}}; makes \code{epsilon} and distances behave +consistently across datasets of different magnitude.} + +\item{targetK}{Integer or NULL. If set, force the clustering result to +exactly \code{targetK} clusters via a seeded K-means refinement (the immune +network supplies the seeds). NULL (default) keeps the emergent cluster count.} + +\item{epsilonQuantile}{Numeric in (0, 1) or NULL. If set, the suppression +threshold adapts each iteration to this quantile of pairwise antibody +distances instead of the fixed \code{epsilon}.} + +\item{verbose}{Logical. If \code{TRUE}, prints progress messages each iteration.} + +\item{...}{Additional arguments forwarded to \code{\link{AINet}}, including +immunology modules (\code{shm}, \code{idiotypic}, \code{germinalCenter}, +\code{microenvironment}, \code{activation}, \code{memory}, +\code{classSwitcher}, \code{init}) and \code{consolidate} / +\code{consolidationSteps}.} } \value{ A list: diff --git a/man/honeycombHIVE.Rd b/man/honeycombHIVE.Rd index eda128e..01d8b3a 100644 --- a/man/honeycombHIVE.Rd +++ b/man/honeycombHIVE.Rd @@ -16,6 +16,7 @@ honeycombHIVE( maxIter = 10, collapseMethod = c("centroid", "medoid", "median", "mode"), minClusterSize = NULL, + smallClusterAction = c("merge", "keep", "drop"), distance = "euclidean", verbose = TRUE, refine = FALSE, @@ -56,8 +57,18 @@ layer.} \item{collapseMethod}{One of "centroid","medoid","median","mode".} -\item{minClusterSize}{Minimum cluster size. Smaller clusters can be -merged/discarded if not NULL.} +\item{minClusterSize}{Minimum cluster size. Clusters with fewer members are +handled per \code{smallClusterAction}. \code{NULL} (default) keeps every +cluster regardless of size.} + +\item{smallClusterAction}{One of \code{"merge"} (default), \code{"keep"}, or +\code{"drop"}, controlling clusters below \code{minClusterSize}. +\code{"merge"} reassigns their members to the nearest surviving prototype so +no observation is lost (the rare cluster is absorbed, not deleted). +\code{"keep"} retains the small cluster as its own prototype (use this to +\emph{preserve} rare populations). \code{"drop"} is the legacy behavior that +discards the cluster and its members from subsequent layers. Ignored when +\code{minClusterSize} is \code{NULL}.} \item{distance}{Distance metric for medoid calculation, e.g. "euclidean".}