context("related_nodes") library(ape) ## R-devel has recently changed the default RNG kind: ## \item The default method for generating from a discrete uniform ## distribution (used in \code{sample()}, for instance) has been ## changed. This addresses the fact, pointed out by Ottoboni and ## Stark, that the previous method made \code{sample()} noticeably ## non-uniform on large populations. See \PR{17494} for a ## discussion. The previous method can be requested using ## \code{RNGkind()} if necessary for reproduction of old results. ## Thanks to Duncan Murdoch for contributing the patch and Gabe ## Becker for further assistance. ## and one should be able to reproduce the Debian error using a current ## version of R-devel (assuming that it builds). ## My suggestion is to disable tests affected by this for the time being. ## ----------------------------------------------------------------------- ## ## set.seed(42) ## # sample bifurcating tree ## bi_tree <- ape::rtree(10) ## bi_tree$tip.label <- paste0("t", 1:10) ## nwk is generated by the above source code in R-3.5.2 with ape v5.2 nwk <- paste0("(((((t1:0.9040313873,t2:0.1387101677):0.5603327462,", "t3:0.9888917289):0.4749970816,(((t4:0.3902034671,t5:0.9057381309):", "0.5142117843,t6:0.4469696281):0.0824375581,(t7:0.7375956178,", "t8:0.8110551413):0.83600426):0.9466682326):0.1174873617,", "t9:0.3881082828):0.9782264284,t10:0.6851697294);") bi_tree <- read.tree(text = nwk) # sample non-bifurcating tree multi_tree <- ape::di2multi(bi_tree, tol=0.5) # bifurcating tree with node names named_bi_tree <- bi_tree named_bi_tree$node.label <- paste0("n", 11:19) # non-bifurcating tree with node names named_multi_tree <- multi_tree named_multi_tree$node.label <- paste0("n", 11:16) empty_tbl <- tibble::tibble( parent=integer(0), node=integer(0), branch.length=numeric(0), label=character(0) ) test_that("conversion to table is reversible", { expect_equal(as.phylo(as_tibble(bi_tree)), bi_tree) expect_equal(as.phylo(as_tibble(multi_tree)), multi_tree) expect_equal(as.phylo(as_tibble(named_bi_tree)), named_bi_tree) expect_equal(as.phylo(as_tibble(named_multi_tree)), named_multi_tree) }) test_that("child works for bifurcating trees", { # a leaf has no children ## expect_equal(child(as_tibble(bi_tree), 1), empty_tbl) expect_equal(nrow(child(as_tibble(bi_tree), 1)), 0) # can find node children expect_equal(child(as_tibble(bi_tree), 19)$node, 7:8) # can find root children expect_equal(child(as_tibble(bi_tree), 11)$node, c(10,12)) }) test_that("child works for non-bifurcating trees", { # a leaf has no children ## expect_equal(child(as_tibble(multi_tree), 1), empty_tbl) expect_equal(nrow(child(as_tibble(multi_tree), 1)), 0) # can find node children expect_equal(child(as_tibble(multi_tree), 12)$node, c(3,9,13,14)) # can find root children expect_equal(child(as_tibble(multi_tree), 11)$node, c(10,12)) }) test_that("offspring works on bifurcating trees", { ## expect_equal(offspring(as_tibble(bi_tree), 1), empty_tbl) expect_equal(nrow(offspring(as_tibble(bi_tree), 1)), 0) expect_equal(offspring(as_tibble(bi_tree), 11)$node, (1:19)[-11]) expect_equal(offspring(as_tibble(bi_tree), 17)$node, c(4:6, 18)) }) test_that("offspring works on non-bifurcating trees", { ## expect_equal(offspring(as_tibble(multi_tree), 1), empty_tbl) expect_equal(nrow(offspring(as_tibble(multi_tree), 1)), 0) expect_equal(offspring(as_tibble(multi_tree), 11)$node, (1:16)[-11]) expect_equal(offspring(as_tibble(multi_tree), 14)$node, c(4:8, 15:16)) }) test_that("parent works for bifurcating trees", { ## expect_equal(parent(as_tibble(bi_tree), 11), empty_tbl) expect_equal(nrow(parent(as_tibble(bi_tree), 11)), 0) expect_equal(parent(as_tibble(bi_tree), 1)$node, 15) expect_equal(parent(as_tibble(bi_tree), 17)$node, 16) }) test_that("parent works for non-bifurcating trees", { ## expect_equal(parent(as_tibble(multi_tree), 11), empty_tbl) expect_equal(nrow(parent(as_tibble(multi_tree), 11)), 0) expect_equal(parent(as_tibble(multi_tree), 8)$node, 16) expect_equal(parent(as_tibble(multi_tree), 14)$node, 12) }) test_that("ancestor works for bifurcating trees", { ## expect_equal(ancestor(as_tibble(bi_tree), 11), empty_tbl) expect_equal(nrow(ancestor(as_tibble(bi_tree), 11)), 0) expect_equal(ancestor(as_tibble(bi_tree), 1)$node, 11:15) expect_equal(ancestor(as_tibble(bi_tree), 17)$node, c(11:13, 16)) }) test_that("ancestor works for non-bifurcating trees", { ## expect_equal(ancestor(as_tibble(multi_tree), 11), empty_tbl) expect_equal(nrow(ancestor(as_tibble(multi_tree), 11)), 0) expect_equal(ancestor(as_tibble(multi_tree), 8)$node, c(11,12,14,16)) expect_equal(ancestor(as_tibble(multi_tree), 14)$node, 11:12) }) test_that("MRCA works for bifurcating trees", { ## 11 is the root node ## expect_equal(MRCA(as_tibble(multi_tree), 11, 5), empty_tbl) expect_equal(nrow(MRCA(as_tibble(multi_tree), 11, 5)), 0) expect_equal(MRCA(as_tibble(bi_tree), 5, 7)$node, 16) ## 16 is ancestor of 5 expect_equal(MRCA(as_tibble(bi_tree), 5, 16)$node, 16) }) test_that("MRCA works for non-bifurcating trees", { ## expect_equal(MRCA(as_tibble(multi_tree), 11, 5), empty_tbl) expect_equal(nrow(MRCA(as_tibble(multi_tree), 11, 5)), 0) expect_equal(MRCA(as_tibble(multi_tree), 5, 7)$node, 14) expect_equal(MRCA(as_tibble(multi_tree), 5, 14)$node, 14) }) test_that("sibling works for bifurcating trees", { ## expect_equal(sibling(as_tibble(bi_tree), 11), empty_tbl) expect_equal(nrow(sibling(as_tibble(bi_tree), 11)), 0) expect_equal(sibling(as_tibble(bi_tree), 1)$node, 2) expect_equal(sibling(as_tibble(bi_tree), 17)$node, 19) }) test_that("sibling works for non-bifurcating trees", { ## expect_equal(sibling(as_tibble(multi_tree), 11), empty_tbl) expect_equal(nrow(sibling(as_tibble(multi_tree), 11)), 0) expect_equal(sibling(as_tibble(multi_tree), 12)$node, 10) expect_equal(sibling(as_tibble(multi_tree), 3)$node, c(9,13,14)) expect_equal(sibling(as_tibble(multi_tree), 4)$node, 5) })