test vs IQTree

data-raw/quartet-concordance-tests.R

@@ -0,0 +1,52 @@
+library("TreeTools")
+tree <- ape::read.tree(text = "(a, (b, (c, (d, ((e, f), (g, h))))));")
+mataset <- matrix(c(0, 0, 0, 0, 1, 1, 1, 1,  0,
+                    0, 0, 1, 1, 1, 1, 1, 1,  0,
+                    0, 0, 1, 1, 1, 1, 2, 2,  0,
+                    1, 0, 0, 0, 1, 1, 1, 1,  0,
+                    1, 0, 0, 0, 0, 1, 1, 1,  0,
+                    0, 0, 0, 0, 1, 1, 2, 2,  0,
+                    0, 0, 1, 1, 2, 2, 3, 3,  0,
+                    0, 1, 2, 3, 0, 1, 2, 3,  0,
+                    0, 1, 2, 0, 0, 2, 2, 3,  0), 9,
+                  dimnames = list(letters[1:9], NULL))
+dat <- MatrixToPhyDat(mataset[1:8, ])
+iqtreeDir <- stop("Path to folder that contains iqtree2.exe")
+treeFile <- paste0(iqtreeDir, "tree.nwk")
+dataFile <- paste0(iqtreeDir, "data.nex")
+charFile <- paste0(iqtreeDir, "data%i.nex")
+
+write.tree(tree, file = treeFile)
+
+# Expected output per character
+vapply(1:9, function(i) {
+  charIFile <- sprintf(charFile, i)
+  write.nexus.data(dat[, i], file = charIFile, format = "STANDARD")
+  dataLines <- readLines(charIFile)
+  dataLines[[5]] <- sub("MISSING=? GAP=- INTERLEAVE=NO ", "", dataLines[[5]],
+                        fixed = TRUE)
+  writeLines(dataLines, charIFile)
+  system2(paste0(iqtreeDir, "iqtree2.exe"),
+          paste("-t", treeFile, "-s", charIFile,
+                "--scf 100000",
+                "-nt 4"
+          ), stdout = NULL)
+  cf <- read.tree(paste0(treeFile, ".cf.tree"))
+
+  as.numeric(cf[["node.label"]][-(1:2)])
+}, double(5))
+
+# Expected output across entire dataset
+write.nexus.data(dat, file = dataFile, format = "STANDARD")
+dataLines <- readLines(dataFile)
+dataLines[[5]] <- sub("MISSING=? GAP=- INTERLEAVE=NO ", "", dataLines[[5]], fixed = TRUE)
+writeLines(dataLines, dataFile)
+system2(paste0(iqtreeDir, "iqtree2.exe"),
+        paste("-t", treeFile, "-s", dataFile,
+              "--scf 1000000",
+              "-nt 4"
+        ), stdout = NULL)
+cf <- read.tree(paste0(treeFile, ".cf.tree"))
+plot(cf, show.node.label = TRUE, use.edge.length = FALSE)
+
+as.numeric(cf[["node.label"]][-(1:2)])

tests/testthat/test-Concordance.R

@@ -49,6 +49,7 @@ test_that("QuartetConcordance(method = minh)", {
   expect_concordance(8, c(NaN, NaN, 0, NaN, NaN))
   expect_concordance(9, c(NaN, 0, 1 / 2, 0, NaN))
 
+  expect_concordance(6, rep(NaN, 5))
   # Values calculated from summing results above
   expect_equal(unname(QuartetConcordance(tree, dat, method = "minh")), 
                c(5 + 3 + 1,
@@ -60,8 +61,10 @@ test_that("QuartetConcordance(method = minh)", {
                   4 + 3 + 2,
                   12 + 8 + 6 + 2 + 2,
                   6 + 2,
-                  4 + 3),
-               tolerance = 0.1)
+                  4 + 3))
+
+  expect_equal(unname(QuartetConcordance(tree, dat, method = "minh")), 
+               c(56.7, 0, 85.4, 0, 62.5), tolerance = 0.01)
 })
 
 test_that("QuartetConcordance() calculates correct values - weighting", {