--- title: "BIDistances: From distance computation to theory-guided selection" output: rmarkdown::html_vignette vignette: > %\VignetteIndexEntry{BIDistances: From distance computation to theory-guided selection} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r, include = FALSE} knitr::opts_chunk$set( collapse = TRUE, comment = "#>", fig.width = 8, fig.height = 4.5 ) ``` ```{r setup} library(BIDistances) ``` # Why BIDistances? Choosing a distance is part of the scientific model, not merely a technical preprocessing decision. A distance determines which observations are regarded as similar, which neighbourhoods are formed, and which structures a distance-based clustering algorithm can recover. `BIDistances` connects four stages of this workflow: 1. **Compute** distances through one dispatcher. 2. **Classify** their mathematical properties and required conditions. 3. **Compare and select** candidates through their distance distributions when distance-based cluster structures are the objective. 4. **Scale** weighted Minkowski calculations across `parallelDist`, an internal multicore implementation, and optional OpenCL kernels. The package also contains the Gene Ontology-derived TF-IDF distance and specialist functions for probability distributions, nearest-neighbour structures, time series, and toroidal coordinates. # Access to distance methods The following release-specific count uses canonical method names accepted by `DistanceMatrix()` for numerical data. Aliases, the `parallelDist` user-defined `custom` route, direct-only functions, and the configurable mixed-data constructions in `manydist` are not counted. ```{r method-count} MethodCount = data.frame( AccessLayer = c( "Package-defined DistanceMatrix routes", "Predefined parallelDist methods", "Additional philentropy identifiers" ), CanonicalRoutes = c(19, 41, 34), Availability = c( "Core", "Core import", "Suggested package" ), stringsAsFactors = FALSE ) knitr::kable(MethodCount, row.names = FALSE) ``` The core installation therefore exposes `r sum(MethodCount$CanonicalRoutes[1:2])` canonical named numerical-data routes. With `philentropy` installed, the total is `r sum(MethodCount$CanonicalRoutes)`. The 34 additional `philentropy` identifiers are obtained after removing exact overlaps and the `chebyshev` and `squared_euclidean` aliases already resolved by the dispatcher. The package-wide choice is broader than this count because `manydist` supplies mixed-type presets and custom component combinations, while functions such as `TWED_Distance()` and `ToroidalEuclidean_Distance()` are intentionally called directly. # Data and Distance conventions - Data: `[1:n,1:d] numerical matrix with n cases, d variables` - Distance: `[1:d,1:d] numerical matrix of distances` Most methods compare rows. The main exceptions are `DTW_Distance()` and `MSMD_Distance()`, which compare time series stored in columns, and `EndresSchindelin_Distance()`, which compares columns as samples of variables. # One dispatcher for different distance families ```{r unified-interface} data("Hepta") Data = Hepta$Data D_euclidean = DistanceMatrix( Data, method = "euclidean" ) D_minkowski_p3 = DistanceMatrix( Data, method = "minkowski", p = 3, threads = 2 ) D_cosine = DistanceMatrix( Data, method = "Cosine_Distance" ) c( EuclideanRows = nrow(D_euclidean), MinkowskiRows = nrow(D_minkowski_p3), CosineRows = nrow(D_cosine) ) ``` The ordinary lower-case method names are first offered to `parallelDist`. Canonical names ending in `_Distance` select native package routes. When a name is unavailable in `parallelDist`, `DistanceMatrix()` can use `philentropy` if that package is installed. Mixed numerical and categorical data are handled through `method = "manydist"`. # Theory-guided comparison of distance distributions The theoretical motivation is the goal of finding **distance-based cluster structures**: distances inside a partition should be smaller than distances between partitions. Under this goal, a multimodal full distance distribution can indicate different regimes of small and large distances. `DistanceDistributions()` operationalizes this exploratory comparison before clustering. `DistanceDistributionAnalysis()` reports the observed Hartigan dip statistic in `SelectionStatistics$DipStatistic`. Candidate vectors with at most 72,000 pairwise values retain `diptest::dip.test()`. Above that threshold, all analyzable distance columns are processed in one integrated C++ batch with an explicit asymptotic p-value calibration and a workspace reused between columns. The full comparison and density visualization are intentionally shown without being executed during vignette construction because they are computationally intensive and exercise optional plotting internals. Run these chunks interactively to reproduce the result. ```{r distance-selection, message = FALSE, warning = FALSE, eval = FALSE} set.seed(42) Selection = DistanceDistributions( Data = Data, DistanceMethods = c( "euclidean", "manhattan", "minkowski", "chord" ), CosineNonParallel = TRUE, CorrelationDist = TRUE, PlotIt = FALSE, PlotSampleSize = 5000 ) Selection$DistanceChoice Selection$SelectionStatistics[, c( "Distance", "DipStatistic", "DipPValue", "BimodalityAmplitude", "Status" )] ``` When `minkowski` is included, the current implementation examines several exponents in addition to the requested base method. Some candidates with `0 < p < 1` are fractional dissimilarities rather than metrics. Their inclusion can be useful for exploration, but their mathematical status must be kept separate from that of Minkowski metrics with `p >= 1`. ```{r distance-selection-plot, message = FALSE, warning = FALSE, eval = FALSE} Selection$ggobject ``` The returned object contains: - `DistanceChoice`: the leading candidate or `"Not Found"`; - `DistanceMatrix`: the selected matrix when a candidate was found; - `OrderedDistances`: the ordered distance features used for comparison; - `ggobject`: the mirrored-density visualization. This procedure addresses a specific modeling goal. Multimodality is evidence for distance-based structure under that goal; it is not a universal proof of cluster validity, biological meaning, or superiority for every downstream algorithm. # Mathematical properties are explicit The word *distance* is used broadly in software. Some methods are metrics, others are pseudometrics, squared metric forms, non-metric dissimilarities, or divergences. `DistanceProperties()` makes these distinctions inspectable. ```{r properties} Properties = DistanceProperties() Relevant = grepl( "Minkowski|Fractional|Tfidf|DTW|Mahalanobis|Cosine", Properties$Function ) knitr::kable( Properties[Relevant, c("Function", "Classification", "Conditions")], row.names = FALSE ) ``` This is important when a downstream algorithm assumes the triangle inequality or identity of indiscernibles. For example: - `Minkowski_Distance()` is a metric for `p >= 1` with strictly positive weights and may become a pseudometric when zero weights collapse distinct vectors. - `Fractional_Distance()` with `0 < p < 1` generally violates the triangle inequality. - `Mahalanobis_Distance()` returns the squared form; applying `sqrt()` yields the metric. - `Cosine_Distance()` returns `1 - cosine similarity` and is not a metric in general. - `Tfidf_Distance()` is a pseudometric on the original annotation rows because distinct rows can have the same scalar TF-IDF weight. # Weighted Minkowski distances on CPU and OpenCL For finite `p >= 1`, the implemented weighted Minkowski distance is \[ d_{p,w}(x,y) = \left(\sum_{k=1}^{d} w_k |x_k-y_k|^p\right)^{1/p}. \] The three backends have different practical roles: | Backend | Role | |---|---| | `parallelDist` | Established multithreaded baseline; weights are represented by coordinate scaling | | `multicore` | Internal `RcppParallel` implementation and CPU fallback | | `opencl` | Optional OpenCL implementation with full, batched-output, and blockwise execution | | `auto` | Tries OpenCL and falls back to the internal multicore implementation with a warning | The Euclidean entry points remain explicit rather than being hidden by the new generalization. `DistanceMatrix(..., method = "euclidean")` uses `parallelDist`, `EuclideanMulticore_Distance()` exposes the internal unweighted CPU implementation, `EuclideanGPU_Distance()` provides the established weighted CPU/OpenCL route, and `ToroidalEuclidean_Distance()` represents periodic two-dimensional geometry. `Minkowski_Distance(..., p = 2)` reuses the optimized Euclidean OpenCL implementation instead of replacing it with a generic power kernel. ## Cross-backend agreement on the CPU ```{r minkowski-cpu} MinkowskiData = as.matrix(iris[1:40, 1:4]) Weights = c(1, 2, 0.5, 1) D_parallel = Minkowski_Distance( Data = MinkowskiData, p = 3, Weights = Weights, backend = "parallelDist", threads = 2 ) D_multicore = Minkowski_Distance( Data = MinkowskiData, p = 3, Weights = Weights, backend = "multicore", threads = 2 ) c( MaximumAbsoluteDifference = max(abs(D_parallel - D_multicore)), EqualWithinTolerance = isTRUE( all.equal(D_parallel, D_multicore, tolerance = 1e-10) ) ) ``` The same route is available through the package dispatcher: ```{r minkowski-dispatch} D_dispatch = DistanceMatrix( MinkowskiData, method = "Minkowski_Distance", p = 3, Weights = Weights, backend = "multicore", threads = 2 ) isTRUE(all.equal(D_dispatch, D_multicore, tolerance = 1e-10)) ``` ## OpenCL route ```{r minkowski-opencl, eval = FALSE} D_opencl = Minkowski_Distance( Data = MinkowskiData, p = 3, Weights = Weights, backend = "opencl", Mem = 2 ) max(abs(D_opencl - D_multicore)) ``` For `p = 2`, the OpenCL route delegates to the established `EuclideanGPU_Distance()` implementation. General exponents use separate Minkowski kernels, so the optimized Euclidean code remains isolated. For `p = Inf`, unit weights are required and the maximum distance is returned. Before a large OpenCL calculation, the memory plan can be inspected without allocating the full matrix: ```{r memory-plan} MemoryPlan = calculateMemoryDemandGPU( n = 70000, d = 784, mem = 6 ) MemoryPlan ``` The GPU backend still returns a dense pairwise matrix, whose storage grows quadratically with the number of cases. The memory planner chooses among complete, output-batched, and input-blocked execution, but it does not remove the quadratic size of the final result. # Gene Ontology-derived TF-IDF distances `Tfidf_Distance()` converts a gene-by-GO-term feature matrix into one scalar TF-IDF weight per gene and then computes absolute pairwise differences between those weights. ```{r tfidf} data("Hearingloss_N109") Gene2Term = Hearingloss_N109$FeatureMatrix_Gene2Term GO_Result = Tfidf_Distance( Gene2Term, tf_fun = mean ) c( Genes = nrow(Gene2Term), GOTerms = ncol(Gene2Term), DistanceRows = nrow(GO_Result$Distance) ) GO_Result$TfidfWeights[1:5] GO_Result$Distance[1:5, 1:5] ``` The distance-only result can be obtained through `DistanceMatrix()`: ```{r tfidf-dispatch} GO_Distance = DistanceMatrix( Gene2Term, method = "Tfidf_Distance", tf_fun = mean ) isTRUE(all.equal(GO_Distance, GO_Result$Distance)) ``` This route illustrates a central advantage of the package: a domain-derived distance, its mathematical classification, and the general distance-distribution workflow are available in one software environment. # Time-series and mixed-data routes `MSMD_Distance()` and `DTW_Distance()` interpret each column as one time series and use row order as temporal order. No time variable should be supplied. ```{r time-series} TimeSeries = cbind( SeriesA = c(0, 1, 2, 1, 0), SeriesB = c(0, 1, 1.5, 1, 0), SeriesC = c(2, 1, 0, 1, 2) ) D_msmd = DistanceMatrix( TimeSeries, method = "MSMD_Distance", ParameterC = 1 ) D_msmd ``` `TWED_Distance()` requires explicit timestamps and is intentionally called directly. This avoids guessing a time variable from the standard matrix interface. When `manydist` is installed, mixed numerical and categorical data can be passed through the same dispatcher: ```{r mixed-data, eval = FALSE} MixedData = data.frame( expression = c(2.1, 1.7, 4.2, 3.8), age = c(42, 38, 61, 58), subtype = factor(c("A", "A", "B", "B")) ) D_mixed = DistanceMatrix( MixedData, method = "manydist", preset = "gower" ) ``` # Which entry point should be used? | Task | Recommended entry point | |---|---| | Compute one of many pairwise matrices | `DistanceMatrix()` | | Compare candidate distance distributions for clustering | `DistanceDistributions()` | | Inspect metric or divergence status | `DistanceProperties()` | | Compute weighted Minkowski distances with selectable backends | `Minkowski_Distance()` | | Compute the optimized weighted Euclidean CPU/OpenCL route | `EuclideanGPU_Distance()` | | Obtain both GO-derived distances and TF-IDF weights | `Tfidf_Distance()` | | Use explicit timestamps for TWED | `TWED_Distance()` | | Compute toroidal coordinates | `ToroidalEuclidean_Distance()` | # References [Thrun, 2021] Thrun, M. C.: The Exploitation of Distance Distributions for Clustering, *International Journal of Computational Intelligence and Applications*, Vol. 20(3), 2150016, DOI: 10.1142/S1469026821500164, 2021. [Stier & Thrun, 2023] Stier, Q., & Thrun, M. C.: Deriving Homogeneous Subsets from Gene Sets by Exploiting the Gene Ontology, *Informatica*, Vol. 34(2), pp. 357-386, DOI: 10.15388/23-INFOR517, 2023. [Aggarwal et al., 2001] Aggarwal, C. C., Hinneburg, A., & Keim, D. A.: On the Surprising Behavior of Distance Metrics in High Dimensional Space, *Database Theory — ICDT 2001*, pp. 420-434, 2001.