Updated data and tests

R/make.demo.data.R

@@ -6,10 +6,11 @@
 #     file.remove("../data/BeckLee_mat50.rda")
 #     file.remove("../data/BeckLee_mat99.rda")
 
+#     set.seed(1)
 #     library(dispRity)
 #     library(paleotree)
 #     library(geiger)
-#     source("../tests/testthat/make.data/multi.ace.R")
+#     source("../tests/testthat/make.data/multi.ace_internal.R")
 #     source("../tests/testthat/make.data/convert.tokens.R")
 #     source("../tests/testthat/make.data/read.nexus.data.R")
 #     ## matrix
@@ -34,9 +35,8 @@
 #     FADLAD <- read.csv("../tests/testthat/make.data/Beck2014_FADLAD.csv", row.names = 1)
 #     FADLAD <- FADLAD[-which(is.na(match(rownames(FADLAD), tree$tip.label))),]
 
-
 #     ## Add the ancestral states estimates
-#     ancestral_states <- multi.ace(matrix, tree, models = "ER", verbose = TRUE)[[1]]
+#     ancestral_states <- multi.ace_internal(matrix, tree, models = "ER", verbose = TRUE)[[1]]
 #     rownames(ancestral_states) <- tree$node.labels
 
 #     ## Combine both
@@ -101,4 +101,8 @@
 
 #     ## save the data
 #     save(BeckLee_disparity, file = "../data/BeckLee_disparity.rda")
-# }
+# }
+
+# # make.demo.data_BeckLeeXXX()
+# # make.demo.data_disparity()
+# # make.demo.data_BeckLee_disparity()

tests/testthat/make.data/make.test.data.R

@@ -1,7 +1,7 @@
 library(dispRity)
 library(paleotree)
 library(geiger)
-source("multi.ace.R")
+source("multi.ace_internal.R")
 source("convert.tokens.R")
 source("read.nexus.data.R")
 
@@ -35,7 +35,9 @@ divRate <- srRes[[1]][1]
 tree <- paleotree::cal3TimePaleoPhy(cladogram, rangesCont, brRate = divRate, extRate = divRate, sampRate = sRate, ntrees = 2, plot = FALSE)
 tree[[1]]$node.label <- tree[[2]]$node.label <- paste0("n", 1:Nnode(tree[[1]]))
 ## Scale the trees to have the same most recent root age
-tree[[1]]$root.time <- tree[[2]]$root.time <- tree[[2]]$root.time
+## Add extra branch length to the root edge
+tree[[1]]$edge.length[which(tree[[1]]$edge[,1] == Ntip(tree[[1]])+1)] <- tree[[1]]$edge.length[which(tree[[1]]$edge[,1] == Ntip(tree[[1]])+1)] + abs(tree[[1]]$root.time - tree[[2]]$root.time)
+tree[[1]]$root.time <- tree[[2]]$root.time
 ## Make the dummy data
 set.seed(1)
 data <- matrix(rnorm((Ntip(tree[[1]])+Nnode(tree[[1]]))*6), nrow = Ntip(tree[[1]])+Nnode(tree[[1]]), ncol = 6, dimnames = list(c(tree[[1]]$tip.label, tree[[1]]$node.label)))
@@ -62,7 +64,7 @@ set.seed(1)
 ## Matches the trees and the matrices
 ## A bunch of trees
 make.tree <- function(n, fun = rtree) {
-    ## Make the tree
+    ## Make the tree
     tree <- fun(n)
     tree <- chronos(tree, quiet = TRUE,
                     calibration = makeChronosCalib(tree, age.min = 10, age.max = 10))
@@ -91,7 +93,7 @@ do.ace <- function(tree, matrix) {
     return(rbind(matrix, apply(matrix, 2, fun.ace, tree = tree)))
 }
 
-## All matrices
+## All matrices
 matrices <- lapply(trees, do.ace, matrix_base)
 
 bound_test_data <- list("matrices" = matrices, "trees" = trees)

tests/testthat/make.data/multi.ace.R → tests/testthat/make.data/multi.ace_internal.R

@@ -54,7 +54,7 @@
 # test <- multi.ace(matrix, tree, models = "ER", use.poly = TRUE, use.uncertain = TRUE, verbose = TRUE)
 
 ##TODO: allow tree to be a multiPhylo object + a sample element that randomly samples a tree everytime and runs ACE on all trees?
-multi.ace <- function(matrix, tree, models, use.poly = FALSE, use.uncertain = FALSE, use.inapp = FALSE, threshold = TRUE, verbose, parallel = FALSE, special.tokens) {
+multi.ace_internal <- function(matrix, tree, models, use.poly = FALSE, use.uncertain = FALSE, use.inapp = FALSE, threshold = TRUE, verbose, parallel = FALSE, special.tokens) {
 
     ## SANITIZING
 
@@ -78,7 +78,7 @@ multi.ace <- function(matrix, tree, models, use.poly = FALSE, use.uncertain = FA
 
     ## Threshold
     #check.class(threshold, c("logical", "numeric"))
-    if(class(threshold) == "logical") {
+    if(is(threshold, "logical")) {
         if(threshold) {
             ## Use the relative threshold function
             threshold.type <- "relative"
@@ -92,16 +92,16 @@ multi.ace <- function(matrix, tree, models, use.poly = FALSE, use.uncertain = FA
     }
 
     #check.class(tree, c("phylo", "multiPhylo"))
-    if(class(tree) == "phylo") {
+    if(is(tree, "phylo")) {
         tree <- list(tree)
         class(tree) <- "multiPhylo"
     }
 
 
     #check.class(matrix, c("matrix", "list"))
     ## Convert the matrix if not a list
-    class_matrix <- class(matrix)
-    if(class_matrix == "list") {
+    class_matrix <- class(matrix)[[1]]
+    if(is(matrix, "list")) {
         matrix <- do.call(rbind, matrix)
     }
 

tests/testthat/test-boot.matrix.R

@@ -458,11 +458,11 @@ test_that("boot.matrix works with multiple trees AND probabilities", {
     test <- boot.matrix(time_slices_multree_proba, bootstraps = 7)
     expect_is(test, "dispRity")
     expect_equal(sort(unlist(lapply(test$subsets, lapply, length), use.name = FALSE)),
-                 c(18, 21, 42, 49, 60, 70))
+                 c(18, 21, 36, 42, 60, 70))
     test <- boot.matrix(time_slices_multree_proba, bootstraps = 7, rarefaction = TRUE)
     expect_is(test, "dispRity")
     expect_equal(sort(unlist(lapply(test$subsets, lapply, length), use.name = FALSE)),
-                 c(18, 21, 21, 21, 28, 28, 35, 35, 42, 42, 42, 49, 49, 56, 60, 63, 70))
+                 c(18, 21, 21, 21, 28, 28, 35, 35, 36, 42, 42, 49, 56, 60, 63, 70))
 
     warn <- capture_warning(boot.matrix(time_slices_multree_proba, bootstraps = 7, boot.type = "single"))
     expect_equal(warn[[1]], "Multiple trees where used in time_slices_multree_proba. The 'boot.type' option is set to \"full\".")

tests/testthat/test-chrono.subsets.R

@@ -572,20 +572,20 @@ test_that("chrono.subsets works with multiPhylo", {
 
     expect_is(test, "dispRity")
     expect_equal(names(test), c("matrix", "tree", "call", "subsets"))
-    expect_equal(names(test$subsets), c("9.31", "4.66", "0"))
+    expect_equal(names(test$subsets), c("9.74", "4.87", "0"))
     expect_equal(unique(unlist(lapply(test$subsets, names), use.names = FALSE)), "elements")
     expect_equal(unlist(lapply(test$subsets, lapply, dim), use.names = FALSE), c(3, 2, 5, 2, 10, 2))
-    expect_equal(unique(c(test$subsets[[2]]$elements)), c(17, 22, 21, 26, NA, 2, 25, 27))
+    expect_equal(unique(c(test$subsets[[2]]$elements)), c(2, 17, 22, 25, 27, 21, 24))
 
     ## Works with discrete
     test <- chrono.subsets(data, tree, method = "discrete", time = 3, inc.nodes = TRUE)
     expect_is(test, "dispRity")
-    expect_equal(unlist(lapply(test$subsets, lapply, dim), use.names = FALSE), c(6, 2, 5, 2, 14, 2))
+    expect_equal(unlist(lapply(test$subsets, lapply, dim), use.names = FALSE), c(5, 2, 7, 2, 14, 2))
 
     ## Works with probabilities
     test <- chrono.subsets(data, tree, method = "continuous", time = 3, model = "gradual.split")
     expect_is(test, "dispRity")
-    expect_equal(unlist(lapply(test$subsets, lapply, dim), use.names = FALSE), c(3, 6, 7, 6, 10, 6))
+    expect_equal(unlist(lapply(test$subsets, lapply, dim), use.names = FALSE), c(3, 6, 6, 6, 10, 6))
 
     ## The output saves the tree
     expect_is(test$tree, "multiPhylo")

tests/testthat/test-dispRity.core.R

@@ -469,9 +469,9 @@ test_that("dispRity deals with probabilities subsets", {
     expect_equal(summary(test2)$n, c(11,20))
     expect_equal(summary(test3)$n, c(15,21))
 
-    expect_equal(as.vector(summary(test1)$obs), c(-0.005, -0.003))
-    expect_equal(as.vector(summary(test2)$obs), c(-0.005, 0.002))
-    expect_equal(as.vector(summary(test3)$obs), c(-0.005, 0.003))
+    expect_equal(as.vector(summary(test1)$obs), c(0.000, 0.002))
+    expect_equal(as.vector(summary(test2)$obs), c(-0.004, 0.000))
+    expect_equal(as.vector(summary(test3)$obs), c(-0.005, 0.002))
 })
 
 test_that("dispRity works with function recycling", {
@@ -540,7 +540,7 @@ test_that("dispRity works with multiple trees from time-slicing", {
     test <- dispRity(boot.matrix(time_slices_proba), metric = c(sum, variances))
     expect_is(test, "dispRity")
     sum_test3 <- summary(test)
-    expect_equal(sum_test3$n, c(3, 7, 10))
+    expect_equal(sum_test3$n, c(3, 6, 10))
     # expect_equal_round(sum_test3$obs.median[c(1,3)], sum_test1$obs.median[c(1,3)])
 
     set.seed(1)
@@ -827,7 +827,6 @@ test_that("dispRity works with the tree component", {
     expect_equal(unlist(c(unname(summary(test)))), c("1:2", "5", "5", "3", "-0.8", "1", "5", "5.9"))
 })
 
-
 test_that("dispRity works with dist.data", {
 
     set.seed(1)

tests/testthat/test-dispRity.metric.R

@@ -353,7 +353,7 @@ test_that("ancestral.dist", {
     data(BeckLee_tree)
     data <- chrono.subsets(BeckLee_mat99, BeckLee_tree, method = "continuous", model = "acctran", time = 5)
     test <- dispRity(data, metric = ancestral.dist)
-    expect_equal(summary(test)$obs.median, c(2.457, 2.538, 2.677, 2.746, 2.741))
+    expect_equal(summary(test)$obs.median, c(2.401, 2.486, 2.621, 2.701, 2.697))
 })
 
 test_that("span.tree.length", {
@@ -719,7 +719,7 @@ test_that("point.dist", {
     data(BeckLee_mat99)
     test <- chrono.subsets(BeckLee_mat99, BeckLee_tree, method = "continuous", model = "equal.split", time = 10)
     test2 <- dispRity(test, metric = point.dist, between.groups = TRUE)
-    expect_equal(summary(test2)$obs.median, c(1.594, 1.838, 1.843, 1.969, 1.828, 1.977, 1.934, 1.892, 1.950))
+    expect_equal(summary(test2)$obs.median, c(1.558, 1.799, 1.804, 1.930, 1.789, 1.940, 1.896, 1.858, 1.912))
 })
 
 test_that("projections", {
@@ -954,7 +954,6 @@ test_that("roudness works", {
     expect_equal_round(test, 0.1776007)
 })
 
-
 test_that("count.neighbours works", {
     set.seed(1)
     dummy_matrix <- matrix(rnorm(50), 5, 10)

tests/testthat/test-dispRity.utilities.R

@@ -428,7 +428,7 @@ test_that("get.disparity", {
         ,names(data$subsets))
     expect_equal(
         round(test[[5]], digit = 5)
-        ,4.09234)
+        ,3.93353)
 
     test <- get.disparity(data, observed = FALSE)
     expect_is(

tests/testthat/test-dist.helper.R

@@ -165,18 +165,18 @@ test_that("works with bootstraps", {
     test <- dispRity(data = bs_data, metric = pairwise.dist)
     check.class(test, "dispRity")
     expect_equal(dim(summary(test)), c(10, 8))
-    expect_equal(summary(test)$obs.median, c(2.472, 2.537, 2.623, 2.723, 2.750, 2.785, 2.841, 2.867, 2.867, 2.867))
+    expect_equal(summary(test)$obs.median, c(2.416, 2.481, 2.567, 2.668, 2.697, 2.729, 2.788, 2.811, 2.811, 2.811))
 
     test <- dispRity(data = bs_data, metric = pairwise.dist, dist.helper = stats::dist)
     check.class(test, "dispRity")
     expect_equal(dim(summary(test)), c(10, 8))
-    expect_equal(summary(test)$obs.median, c(2.472, 2.537, 2.623, 2.723, 2.750, 2.785, 2.841, 2.867, 2.867, 2.867))
+    expect_equal(summary(test)$obs.median, c(2.416, 2.481, 2.567, 2.668, 2.697, 2.729, 2.788, 2.811, 2.811, 2.811))
 
     dist_matrix <- dist(BeckLee_mat99)
     test <- dispRity(data = bs_data, metric = pairwise.dist, dist.helper = dist_matrix)
     check.class(test, "dispRity")
     expect_equal(dim(summary(test)), c(10, 8))
-    expect_equal(summary(test)$obs.median, c(2.472, 2.537, 2.623, 2.723, 2.750, 2.785, 2.841, 2.867, 2.867, 2.867))
+    expect_equal(summary(test)$obs.median, c(2.416, 2.481, 2.567, 2.668, 2.697, 2.729, 2.788, 2.811, 2.811, 2.811))
 
 
 # test <- microbenchmark("no help"       = dispRity(bs_data, metric = pairwise.dist),

tests/testthat/test-plot.dispRity.R

@@ -40,7 +40,7 @@ test_that("get.plot.params works", {
     ## The plotting options
     expect_equal(plot_params$options$xlab, "Time (Mya)")
     expect_equal(plot_params$options$ylab, "c(median, centroids)")
-    expect_equal_round(plot_params$options$ylim, c(1.546577, 2.012542), 6)
+    expect_equal_round(plot_params$options$ylim, c(1.516207, 1.971640), 6)
     expect_equal(plot_params$options$col, c("black", "#BEBEBE", "#D3D3D3"))
     ## Observed data
     expect_equal(names(plot_params$observed_args), c("observed", "col", "names", "data", "pch", "cex"))

tests/testthat/test-print.dispRity.R

@@ -270,8 +270,8 @@ test_that("print.dispRity with model.test data", {
          "Call: model.test(data = model_test_data, model = models, time.split = 65, fixed.optima = TRUE, verbose = FALSE) ",
          "",
          "        aicc delta_aicc weight_aicc",
-         "BM -31.29071   0.000000   0.7856167",
-         "OU -28.69331   2.597407   0.2143833",
+         "BM -32.69195   0.000000   0.7856167",
+         "OU -30.09454   2.597407   0.2143833",
          "",
          "Use x$full.details for displaying the models details",
          "or summary(x) for summarising them."
@@ -302,7 +302,7 @@ test_that("print.dispRity with model.test data", {
           "",
           "Model simulated (10 times):",
           "    aicc log.lik param ancestral state sigma squared",
-          "BM -31.3   17.92     2           3.099         0.002",
+          "BM -32.7   18.62     2           2.967         0.002",
           "",
           "Rank envelope test:",
           " p-value of the global test: 0.1818182 (ties method: erl)",

tests/testthat/test-summary.dispRity.R

@@ -24,12 +24,12 @@ test_that("get.summary", {
     expect_is(test, "list")
     expect_equal(names(test), c("cent_tend", "quantiles"))
     expect_equal(round(test[[1]], digit = 5), round(mean(unlist(disparity$disparity$`30`[[2]])), digit = 5))
-    expect_equal(round(test[[2]], digit = 2), c("25%" = 1.79, "75%" = 1.86))
+    expect_equal(round(test[[2]], digit = 2), c("25%" = 1.75, "75%" = 1.82))
 
     test_no_cent_tend <- get.summary(disparity$disparity$`30`[[2]], quantiles = c(50))
     expect_is(test_no_cent_tend, "list")
     expect_equal(names(test_no_cent_tend), "quantiles")
-    expect_equal(round(test_no_cent_tend[[1]], digit = 2), c("25%" = 1.79, "75%" = 1.86))
+    expect_equal(round(test_no_cent_tend[[1]], digit = 2), c("25%" = 1.75, "75%" = 1.82))
 
     test_no_quant <- get.summary(disparity$disparity$`30`[[2]], cent.tend = mean)
     expect_is(test_no_quant, "list")