update tests for Bruvo

These are based on the values from polysat

src/poppr_distance.c

@@ -313,9 +313,13 @@ all.equal(run_models(tg, "poppr_bruvo"), run_models(tg, "Rbruvo"))
 # works with recursion
 tg1[1, ] <- c(0, 0, 0, 0, 24)
 all.equal(run_models(tg1, "poppr_bruvo"), run_models(tg1, "Rbruvo"))
-all.equal(run_models(tg1, "poppr_bruvo", old_model = TRUE), 
-          run_models(tg1, "Rbruvo", old_model = TRUE))
+all.equal(run_models(tg1, "poppr_bruvo", old_model = TRUE), run_models(tg1, "Rbruvo", old_model = TRUE))
 
+tg2 <- tg1
+tg2[1, -(1:3)] <- c(102, 104)/2
+tg2[2, -1] <- c(104,104,106,110)/2
+all.equal(run_models(tg1, "poppr_bruvo"), run_models(tg1, "Rbruvo"))
+all.equal(run_models(tg1, "poppr_bruvo", old_model = TRUE), run_models(tg1, "Rbruvo", old_model = TRUE))
 
 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
 SEXP bruvo_distance(SEXP bruvo_mat, SEXP permutations, SEXP alleles, SEXP m_add, SEXP m_loss, SEXP old_model)
@@ -466,29 +470,6 @@ int fact(int x)
 	woo: p * p!
 	loss: TRUE/FALSE: impute under genome loss model.
 	add: TRUE/FALSE: impute under genome addition model. 
-
-	Test code comparing current status to polysat's Bruvo2.distance:
-================================================================================
-poppr_bruvo <- function(){ 
-  return(c(.Call("single_bruvo", c(20,23,24,0,20,24,26,43), .Call("permuto", 4), 4, 0, 0),
-.Call("single_bruvo", c(20,23,24,0,20,24,26,43), .Call("permuto", 4), 4, 1, 0),
-.Call("single_bruvo", c(20,23,24,0,20,24,26,43), .Call("permuto", 4), 4, 0, 1),
-.Call("single_bruvo", c(20,23,24,0,20,24,26,43), .Call("permuto", 4), 4, 1, 1)
-))
-}
-
-polysat_bruvo <- function(){
-  return(c(Bruvo2.distance(c(20,23,24), c(20,24,26,43), usatnt=1, loss=FALSE, add=FALSE),
-Bruvo2.distance(c(20,23,24), c(20,24,26,43), usatnt=1, loss=FALSE, add=TRUE),
-Bruvo2.distance(c(20,23,24), c(20,24,26,43), usatnt=1, loss=TRUE, add=FALSE),
-Bruvo2.distance(c(20,23,24), c(20,24,26,43), usatnt=1, loss=TRUE, add=TRUE)
-))
-}
-
-library(polysat)
-polysat_bruvo()
-poppr_bruvo()
-polysat_bruvo() == poppr_bruvo()
 ==============================================================================*/
 double bruvo_dist(int *in, int *nall, int *perm, int *woo, int *loss, int *add, int old_model)
 {

tests/testthat/test-values.R

@@ -89,9 +89,22 @@ test_that("Bruvo's distance will go through the recusion", {
   addLOSS <- as.vector(bruvo.dist(testgid, add = FALSE, loss = TRUE))
   ADDLOSS <- as.vector(bruvo.dist(testgid, add = TRUE, loss = TRUE))
   expect_equal(ADDloss, 0.671874523162842)
-  expect_equal(addLOSS, 0.351561918854713)
+  expect_equal(addLOSS, 0.293456600047648)
   expect_equal(addloss, 0.75)
-  expect_equal(ADDLOSS, 0.511718221008778)
+  expect_equal(ADDLOSS, 0.482665561605245)
+})
+test_that("Multinomial coefficient respects index, not value", {
+  skip_on_cran()
+  testdf  <- data.frame(test = c("00/00/00/51/52", "00/52/52/53/55"))
+  testgid <- df2genind(testdf, ploidy = 5, sep = "/")
+  ADDloss <- as.vector(bruvo.dist(testgid, add = TRUE, loss = FALSE))
+  addloss <- as.vector(bruvo.dist(testgid, add = FALSE, loss = FALSE))
+  addLOSS <- as.vector(bruvo.dist(testgid, add = FALSE, loss = TRUE))
+  ADDLOSS <- as.vector(bruvo.dist(testgid, add = TRUE, loss = TRUE))
+  expect_equal(ADDloss, 0.4375)
+  expect_equal(addLOSS, 0.25)
+  expect_equal(addloss, 0.625)
+  expect_equal(ADDLOSS, 0.34375)
 })
 
 test_that("Repeat lengths can be in any order and length if named", {