suppressPackageStartupMessages(library(tidyverse))
md <- read.csv("../../data/Haws/sample_metadata.csv") %>% dplyr::select(-1) %>% arrange(sample_number) %>% rename(pH = ph) %>% mutate(cell = paste(ploidy, pH, sep = "-"))
ids <- c(paste0("2H-",1:6),paste0("2L-",1:6),paste0("3H-",1:6),paste0("3L-",1:6))
perm1 <- function(metadata, perm) { set.seed(perm)
  metadata %>% group_by(pH) %>% mutate(ploidy = sample(ploidy)) %>% ungroup() %>% group_by(ploidy) %>% mutate(pH = sample(pH)) %>% ungroup() }
counts <- read_csv("permutation-counts-cov5-hiperc99.9-cell5-meth10-All-interaction.csv", show_col_types = FALSE) %>% filter(fdr == 0.05, minDiff == 0) %>% dplyr::select(perm, term, DML) %>% pivot_wider(names_from = term, values_from = DML)
out <- map_dfr(1:20, function(p) { d <- perm1(md, p)
  tibble(perm = p, ploidyAgree = sum(d$ploidy == md$ploidy), pHAgree = sum(d$pH == md$pH),
         g3H2 = paste(d$ploidy[14], d$pH[14]), pc2 = paste(paste0(substr(d$ploidy[c(1,2,17)],1,1), substr(d$pH[c(1,2,17)],1,1)), collapse = " ")) })
res <- left_join(out, counts, by = "perm")
print(as.data.frame(res))
cat("cor(|ploidyAgree-12|, ploidy DML):", cor(abs(res$ploidyAgree - 12), res$ploidy, method = "spearman"), "\n")
cat("cor(|pHAgree-12|, pH DML):", cor(abs(res$pHAgree - 12), res$pH, method = "spearman"), "\n")