test-bitwise.R
In poppr: Genetic Analysis of Populations with Mixed Reproduction

context("Bitwise distance calculations")
set.seed(991)
dat <- sample(c(0, 1, NA), 50, replace = TRUE, prob = c(0.49, 0.5, 0.01))
mat <- matrix(dat, nrow = 5, ncol = 10)
mat2 <- mat
mat2[-1, ] <- mat2[-1, ] * 2

test_that("bitwise.dist works for haploids", {
  mat.gl <- new("genlight", mat, parallel = FALSE)
  expected <- rowSums(apply(mat, 2, dist), na.rm = TRUE)/10
  expect_equivalent(as.vector(bitwise.dist(mat.gl, threads = 1L)), expected)
  expect_true(all(ploidy(mat.gl) == 1))
})

test_that("bitwise.dist works for diploids with only one type of homozygote", {
  mat2.gl <- new("genlight", mat2, parallel = FALSE)
  expect_error(bitwise.dist(mat2.gl, threads = 1L), "ploidy")
  ploidy(mat2.gl) <- rep(2, 5)
  expected <- rowSums(apply(mat2, 2, dist), na.rm = TRUE)/20
  expect_equivalent(as.vector(bitwise.dist(mat2.gl, threads = 1L)), expected)
})

test_that("bitwise.dist can do euclidean", {
  mat2.gl <- new("genlight", mat2, parallel = FALSE)
  ploidy(mat2.gl) <- rep(2, 5)
  expect_equivalent(bitwise.dist(mat2.gl, scale_missing = TRUE, euclid = TRUE, threads = 1L), dist(mat2.gl))
})

test_that("bitwise.dist can do euclidean with lots of missing data", {
  # skip_on_cran()
  set.seed(999)
  mat2[sample(length(mat2), 10)] <- NA
  mat2.gl <- new("genlight", mat2, parallel = FALSE)
  ploidy(mat2.gl) <- rep(2, 5)
  expect_equivalent(bitwise.dist(mat2.gl, scale_missing = TRUE, euclid = TRUE, threads = 1L), dist(mat2.gl))
})

test_that("bitwise.dist can actually handle genind objects", {
  # skip_on_cran()
  data("partial_clone", package = "poppr")
  pdist <- diss.dist(partial_clone, percent = TRUE)
  expect_equivalent(pdist, bitwise.dist(partial_clone))
})


test_that("bitwise.dist produces reasonable results for diploids", {

  # skip_on_cran()
  dat <- list(c(2,2,2,2,2,2,2,2,2,0),
              c(1,1,1,0,0,0,0,0,0,2),
              c(2,2,2,2,2,2,2,2,2,2),
              c(2,2,2,2,2,2,2,2,2,0),
              c(2,NA,NA,NA,NA,NA,NA,NA,NA,NA))
  z <- new("genlight",dat, parallel = FALSE)
  
  # Sun May 13 12:37:19 2018 ------------------------------
  # Jonah originally set this up so it was simple to test, but did not account
  # for the fact that we needed to test when the missing data were in both the
  # i and the j positions. In this case, the missing are only in the j position.
  # Luckily, we can easily test this by inverting the data and results. In order
  # to distinguish between the two forms, I'll append i and j.
  missing_match_dif_i <- bitwise.dist(z[5:1, ], missing_match = TRUE, mat = TRUE, differences_only = TRUE)
  missing_match_dif_j <- bitwise.dist(z, missing_match = TRUE, mat = TRUE, differences_only = TRUE)
  
  expected_match_dif <- c(0.0, 1.0, 0.1, 0.0, 0.0, 
                          1.0, 0.0, 0.9, 1.0, 0.1, 
                          0.1, 0.9, 0.0, 0.1, 0.0, 
                          0.0, 1.0, 0.1, 0.0, 0.0, 
                          0.0, 0.1, 0.0, 0.0, 0.0)
  dim(expected_match_dif) <- c(5L, 5L)
  
  missing_nomatch_dif_i <- bitwise.dist(z[5:1, ], missing_match = FALSE, mat = TRUE, differences_only = TRUE)
  missing_nomatch_dif_j <- bitwise.dist(z, missing_match = FALSE, mat = TRUE, differences_only = TRUE)
  
  expected_nomatch_dif <- c(0.0, 1.0, 0.1, 0.0, 0.9, 
                            1.0, 0.0, 0.9, 1.0, 1.0, 
                            0.1, 0.9, 0.0, 0.1, 0.9, 
                            0.0, 1.0, 0.1, 0.0, 0.9, 
                            0.9, 1.0, 0.9, 0.9, 0.0)
  dim(expected_nomatch_dif) <- c(5L, 5L)
  
  expect_equivalent(missing_match_dif_i,   expected_match_dif[5:1, 5:1])
  expect_equivalent(missing_nomatch_dif_i, expected_nomatch_dif[5:1, 5:1])
  
  expect_equivalent(missing_match_dif_j,   expected_match_dif)
  expect_equivalent(missing_nomatch_dif_j, expected_nomatch_dif)
  
  missing_match_dist_i <- bitwise.dist(z[5:1, ], missing_match = TRUE, mat = TRUE, differences_only = FALSE)
  missing_match_dist_j <- bitwise.dist(z, missing_match = TRUE, mat = TRUE, differences_only = FALSE)
  expected_match_dist <- c(0.0, 17.0/20.0, 2.0/20.0, 0.0, 0.0, 
                           17.0/20.0, 0.0, 15.0/20.0, 17.0/20.0, 1.0/20.0, 
                           2.0/20.0, 15.0/20.0, 0.0, 2.0/20.0, 0.0, 
                           0.0, 17.0/20.0, 2.0/20.0, 0.0, 0.0, 
                           0.0, 1.0/20.0, 0.0, 0.0, 0.0)
  dim(expected_match_dist) <- c(5L, 5L)
  expect_equivalent(missing_match_dist_i, expected_match_dist[5:1, 5:1])
  expect_equivalent(missing_match_dist_j, expected_match_dist)
})

test_that("bitwise.ia produce reasonable results for haploids", {
  # skip_on_cran()
  dat <- list(c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1),
              c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
              c(1, 1, NA, NA, NA, NA, NA, NA, NA, NA))
  z   <- new("genlight",dat, parallel = FALSE)
  expected_dist      <- c(0, 1, 0, 1, 0, 0.2, 0, 0.2, 0)
  dim(expected_dist) <- c(3L, 3L)
  
  missing_match_dif_i <- bitwise.dist(z[3:1, ], missing_match = TRUE, mat = TRUE, differences_only = TRUE)
  missing_match_dif_j <- bitwise.dist(z, missing_match = TRUE, mat = TRUE, differences_only = TRUE)
  
  expect_equivalent(missing_match_dif_i, expected_dist[3:1, 3:1])
  expect_equivalent(missing_match_dif_j, expected_dist)
  
  expected_nomatch_dist <- c(0, 1, .8, 1, 0, 1, .8, 1, 0)
  dim(expected_nomatch_dist) <- c(3L, 3L)

  missing_nomatch_dif_i <- bitwise.dist(z[3:1, ], missing_match = FALSE, mat = TRUE, differences_only = TRUE)
  missing_nomatch_dif_j <- bitwise.dist(z, missing_match = FALSE, mat = TRUE, differences_only = TRUE)
  
  expect_equivalent(missing_nomatch_dif_i, expected_nomatch_dist[3:1, 3:1])
  expect_equivalent(missing_nomatch_dif_j, expected_nomatch_dist)
  
})

context("bitwise.ia cromulence")

test_that("bitwise.ia can use both missing-match and missing-nomatch ", {
  # skip_on_cran()
  dat <- list(c(2,2,2,2,2,2,2,2,2,0),
              c(1,1,1,0,0,0,0,0,0,2),
              c(2,2,2,2,2,2,2,2,2,2),
              c(2,2,2,2,2,2,2,2,2,0),
              c(2,NA,NA,NA,NA,NA,NA,NA,NA,NA))
  z <- new("genlight",dat, parallel = FALSE)
  tz <- tab(z, NA.method = "asis")
  tz[tz == 2] <- "A/A"
  tz[tz == "1"] <- "A/B"
  tz[tz == "0"] <- "B/B"
  tzg <- df2genind(tz, sep = "/")
  
  # Missing-match is the default and matches that of the non-bitwise version
  expect_equal(bitwise.ia(z, missing_match = TRUE), ia(tzg)[[2]])
  expect_equal(bitwise.ia(z[5:1, ], missing_match = TRUE), ia(tzg)[[2]])
  
  ianm <- function(i, z){
    as.vector(bitwise.dist(z[, i], percent = FALSE, missing_match = FALSE))
  }
  np  <- choose(5, 2)
  bdl <- vapply(seq(10), ianm, integer(np), z)
  dlist <- list(
             d.vector = colSums(bdl),
             d2.vector = colSums(bdl * bdl),
             D.vector = rowSums(bdl)
           )
  res <- poppr:::ia_from_d_and_D(dlist, np)
  expect_equal(res[[2]], bitwise.ia(z, missing_match = FALSE))
})