Starting up

We will need to load some libraries:

library(ape)
library(wiritttes)
library(ggplot2)
library(nLTT)
library(phangorn)
options(warn = 2)

Creating parameter files

Here, we create the four parameter filenames:

filename <- "add_posterior_demo.RDa"

Now the parameter files are created with the desired parameters:

n_beast_runs <- 2

save_parameters_to_file(
  rng_seed = 2,
  sirg = 0.5,
  siri = 0.5,
  scr = 1.0e-1,
  erg = 0.1,
  eri = 0.1,
  age = 5,
  mutation_rate = 0.01,
  n_alignments = 1,
  sequence_length = 1000,
  nspp = 10,
  n_beast_runs = n_beast_runs,
  filename = filename
)  
testit::assert(is_valid_file(filename))

Running the simulation

The simulation has multiple steps. The last step is adding the posterior to the parameter file.

do_simulation(filename)

Show the posterior

posteriors <- get_posteriors(read_file(filename))
testit::assert(length(posteriors) == 2 * n_beast_runs)
testit::assert(tracerer::is_posterior(wiritttes::get_posterior(file = read_file(filename), sti = 1, ai = 1, pi = 1)))
testit::assert(tracerer::is_posterior(wiritttes::get_posterior(file = read_file(filename), sti = 1, ai = 1, pi = 2)))
testit::assert(tracerer::is_posterior(wiritttes::get_posterior(file = read_file(filename), sti = 2, ai = 1, pi = 1)))
testit::assert(tracerer::is_posterior(wiritttes::get_posterior(file = read_file(filename), sti = 2, ai = 1, pi = 2)))

phylogenies_1 <- wiritttes::get_posterior(file = read_file(filename), sti = 1, ai = 1, pi = 1)$trees
phylogenies_2 <- wiritttes::get_posterior(file = read_file(filename), sti = 1, ai = 1, pi = 2)$trees
phylogenies_3 <- wiritttes::get_posterior(file = read_file(filename), sti = 2, ai = 1, pi = 1)$trees
phylogenies_4 <- wiritttes::get_posterior(file = read_file(filename), sti = 2, ai = 1, pi = 2)$trees

testit::assert(length(phylogenies_1) > 0)
testit::assert(length(phylogenies_2) > 0)
testit::assert(length(phylogenies_3) > 0)
testit::assert(length(phylogenies_4) > 0)

Here another hack: to get the densiTree function working, phylogenies must be of class multiphylo. Let it be so:

class(phylogenies_1) <- "multiPhylo"
class(phylogenies_2) <- "multiPhylo"
class(phylogenies_3) <- "multiPhylo"
class(phylogenies_4) <- "multiPhylo"

Displaying first posterior of the youngest species tree:

phangorn::densiTree(
  phylogenies_1, 
  type = "cladogram", 
  alpha = 1
)

Displaying second posterior of the youngest species tree:

phangorn::densiTree(
  phylogenies_2, 
  type = "cladogram", 
  alpha = 1
)

Displaying first posterior of the oldest species tree:

phangorn::densiTree(
  phylogenies_3, 
  type = "cladogram", 
  alpha = 1
)

Displaying second posterior of the oldest species tree:

phangorn::densiTree(
  phylogenies_4, 
  type = "cladogram", 
  alpha = 1
)
file.remove(filename)
file.remove(dir(pattern="add_posterior_demo_"))


richelbilderbeek/wiritttes documentation built on May 27, 2019, 8:14 a.m.