Docker-in-Docker (DinD) capabilities of public runners deactivated. More info

rer_converge.R 1.91 KB
Newer Older
1 2 3 4 5
suppressPackageStartupMessages({
  library(tidyverse)
  library(RERconverge)
})

6 7
master_tree <- read.tree(file="<<<dep master_tree>>>")
gene_trees <- readTrees("<<<dep gene_tree_set>>>", max.read=<<<max_read>>>, master_tree=master_tree, minSpecs=<<<min_specs>>>)
8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26

relative_rates <- getAllResiduals(
  gene_trees,
  transform = "<<<template_transform transform>>>",
  weighted = <<<template_bool weighted>>>,
  scale = <<<template_bool scale>>>
)

load_species_file <- function(path) {
    read_table(path, col_names = "species") %>% pull(species)
}

load_traits_file <- function(path, is_continuous) {
  if (is_continuous)
    read_table(path, col_names = c('species', 'value')) %>% deframe()
  else
    load_species_file(path)
}

27 28 29 30
phenotypes <- load_traits_file("<<<dep phenotypes>>>", is_continuous)

build_traits_paths <-
  function(phenotypes,
31 32 33 34 35 36
           gene_trees,
           is_continuous,
           useSpecies = NULL,
           clade = c("ancestral", "terminal", "all"))
  {
    if (is_continuous) {
37
      char2Paths(phenotypes, gene_trees)
38
    } else {
39
      foreground2Tree(phenotypes,
40 41 42 43 44 45 46 47 48
                      gene_trees,
                      clade = clade,
                      useSpecies = useSpecies,
                      plot=F # FIXME : plot=T raises cryptic error
                      ) %>%
        tree2Paths(gene_trees)
    }
  }

49 50
traits_paths <- build_traits_paths(
  phenotypes,
51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
  gene_trees,
  <<<template_bool continuous>>>
)

res <- if (<<<template_bool continuous>>>) {
  correlateWithContinuousPhenotype(
    relative_rates,
    traits_paths,
  )
} else {
  correlateWithBinaryPhenotype(
    relative_rates,
    traits_paths,
  )
}

67 68 69 70
res <- res %>% 
  mutate(gene_id=names(gene_trees$trees)) %>% 
  select(gene_id, everything())

71 72 73 74 75 76 77 78
saveRDS(list(
  relative_rates = relative_rates,
  results = res,
  traits_paths = traits_paths,
  master_tree = master_tree,
  gene_trees = gene_trees,
  phenotypes = phenotypes
), file="<<<dest>>>")