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

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

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

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)
}

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

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

50 51
traits_paths <- build_traits_paths(
  phenotypes,
52
  gene_trees,
53
  is_continuous
54 55
)

56
res <- if (is_continuous) {
57 58 59 60 61 62 63 64 65 66 67
  correlateWithContinuousPhenotype(
    relative_rates,
    traits_paths,
  )
} else {
  correlateWithBinaryPhenotype(
    relative_rates,
    traits_paths,
  )
}

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

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