Bugfixes for dist.method='slingshot' choice.

Bugfixes for correct endpoint search.
LTLA · Dec 15, 2020 · 4d8a020 · 4d8a020
1 parent 7307b54
commit 4d8a020
Show file tree

Hide file tree

Showing 3 changed files with 31 additions and 21 deletions.
diff --git a/NAMESPACE b/NAMESPACE
@@ -46,6 +46,7 @@ importFrom(Matrix,crossprod)
 importFrom(Matrix,rowMeans)
 importFrom(Matrix,rowSums)
 importFrom(Matrix,t)
+importFrom(S4Vectors,head)
 importFrom(SingleCellExperiment,colLabels)
 importFrom(SingleCellExperiment,reducedDim)
 importFrom(SummarizedExperiment,assay)

diff --git a/R/createClusterMST.R b/R/createClusterMST.R
@@ -221,7 +221,7 @@ NULL
                 # Distances not really intepretable as Mahalanobis distances anymore.
                 warning("'use.median=TRUE' with 'dist.method=\"", dist.method, "\"' may yield unpredictable results")
             }
-            use.full <- (dist.method == "scaled.full" || (dist.method == "slingshot" && min(table(clusters)) <= ncol(x)))
+            use.full <- (dist.method == "scaled.full" || (dist.method == "slingshot" && min(table(clusters)) > ncol(x)))
             dmat <- .dist_clusters_scaled(x, clusters, centers=centers, full=use.full)
         }
     }
@@ -320,6 +320,7 @@ NULL
     (distances + t(distances))
 }
 
+#' @importFrom S4Vectors head
 .enforce_endpoints <- function(dmat, endpoints, allow.dyads=FALSE) {
     available <- dmat[unique(endpoints),,drop=FALSE]
     best.stats <- new.env()
@@ -335,25 +336,25 @@ NULL
             return(NULL)
         } else if (distance > best.stats$distance) {
             return(NULL)
+        } 
+
+        current <- rownames(available)[i]
+        self.used <- which(path == current)
+
+        if (length(self.used) == 1) { 
+            # Endpoint-to-endpoint dyads should be reciprocated,
+            # with no distance added (if they are allowed to exist).
+            reciprocal <- rownames(available)[self.used]
+            if (!reciprocal %in% path && allow.dyads) {
+                SEARCH(c(path, reciprocal), distance)
+            }
         } else {
-            current <- rownames(available)[i]
-            used <- which(path == current)
-
-            if (length(used) > 1) {
-                # Can't have an endpoint connected to two things.
-                return(NULL)
-            } else if (length(used) == 1) { 
-                # Endpoint-to-endpoint dyads should be reciprocated,
-                # with no distance added.
-                reciprocal <- rownames(available)[used]
-                if (!reciprocal %in% path && allow.dyads) {
-                    SEARCH(c(path, reciprocal), distance)
-                }
-            } else {
-                allowed <- setdiff(colnames(available), c(current, path))
-                for (j in allowed) {
-                    SEARCH(c(path, j), distance + available[i,j])
-                }
+            used.endpoints <- c(current, # currently in use.
+                head(rownames(available), length(path)), # endpoints connected from in previous steps.
+                intersect(path, rownames(available))) # endpoints connected to in previous steps.
+            allowed <- setdiff(colnames(available), used.endpoints)
+            for (j in allowed) {
+                SEARCH(c(path, j), distance + available[i,j])
             }
         }
     }

diff --git a/tests/testthat/test-create-mst.R b/tests/testthat/test-create-mst.R
@@ -135,6 +135,12 @@ test_that("MST construction works with endpoint specification", {
     out <- createClusterMST(y, endpoint="B", clusters=NULL)
     expect_true(igraph::degree(out, "B")==1L)
     expect_true(is.infinite(igraph::E(out)$gain[1]))
+    expect_identical(igraph::components(out)$no, 1L)
+
+    out <- createClusterMST(y, endpoint=c("C", "D"), clusters=NULL)
+    expect_true(igraph::are_adjacent(out, "B", "C")) 
+    expect_true(igraph::are_adjacent(out, "B", "D"))
+    expect_identical(igraph::components(out)$no, 1L)
 
     # Does sensible things when endpoints= doesn't have an effect.
     ref <- createClusterMST(y, clusters=NULL)
@@ -254,6 +260,7 @@ test_that("MST construction works with scaled distances", {
     mst3 <- createClusterMST(y0.simple, clusters=clusters.simple, dist.method="slingshot")
     expect_identical(ref[] > 0, mst3[] > 0)
     expect_identical(igraph::V(ref)$coordinates, igraph::V(mst3)$coordinates)
+    expect_identical(mst2[], mst3[]) # chooses scaled.full.
 
     # Trying something that requires a bit more... finesse.
     ref <- createClusterMST(y.complex, clusters=clusters.complex)
@@ -265,8 +272,9 @@ test_that("MST construction works with scaled distances", {
     mst <- createClusterMST(y.complex, clusters=clusters.complex, dist.method="scaled.full")
     expect_true(igraph::are_adjacent(mst, "1", "2"))
 
-    mst <- createClusterMST(y.complex, clusters=clusters.complex, dist.method="slingshot")
-    expect_true(igraph::are_adjacent(mst, "1", "2"))
+    mst2 <- createClusterMST(y.complex, clusters=clusters.complex, dist.method="slingshot")
+    expect_true(igraph::are_adjacent(mst2, "1", "2"))
+    expect_identical(mst2[], mst[]) # chooses scaled.full.
 
     # Works correctly with weight matrices.
     mat <- factor2matrix(clusters.complex)