TreeProcesses.jl

Some (common) forward and backward processes to generate binary trees as seen for example in population genetics.

Based on BinaryTrees.jl, which in turn implements the interface from AbstractTrees.jl

Installation

The package is unregistered. Install it directly from GitHub.

julia> ]add https://github.com/skleinbo/BinaryTrees.jl https://github.com/skleinbo/TreeProcesses.jl

Currently implemented

Forward:

• birth_death(n, T, d, b=1.0; N=0): Starting from n nodes, a randomly selected one splits with probability b. Independently, a node dies with probability d. Run for T time steps (a birth plus death event are one time step), or until N nodes are present unless N==0 (default). Lineages that die out are automatically pruned, i.e. after sufficiently many (O(n^2)) time steps, a most recent common ancestor will be found. Returns all surviving root nodes.
• maximally_balanced(n)
• maximally_unbalanced(n)
• moran(n,T): birth_death with b/d probabilities both equal to 1 to maintain constant population size.
• yule(n): Starting from a root node, split leaves uniformly n times.

Backward:

• coalescent(n): Neutral coalescent process. Starting from n nodes, pick two at random and merge until only one is left.
• weighted_coalescent(n, w; fuse=max): Like coalescent, but each node is assigned the respective weight from the vector w with which it partakes in a coalescent event. Ancestral nodes' weights are computed from their children's weights by applying fuse. The standard coalescent is recovered by choosing all weights equal and setting fuse=first.

Observables

• AC!(t): Annotate each node of a tree with

  • $A$: Number of nodes in the subtree, including itself.
  • $C$: Cumulative number of nodes in the subtree, i.e. $\sum A(i)$ for $i$ in the subtree, including the node itself.

  Note: The val field of the nodes must be a Vector that can store Ints, which is the default.

Example

julia> using TreeProcesses

julia> T = weighted_coalescent(2^4, rand(2^4))
BinaryTree{Vector{Int64}}([0, 0]) 11038476695863080329 with 2 children and no parent.

# calculate and return number of nodes, A & C
julia> AC!(T)
(31, [31, 29, 27, 23, 21, 7, 3, 1, 1, 3  …  3, 1, 1, 1, 1, 3, 1, 1, 1, 1], [203, 171, 141, 109, 85, 17, 5, 1, 1, 5  …  5, 1, 1, 1, 1, 5, 1, 1, 1, 1])

# the tree has been annotated with the observable values
julia> T
BinaryTree{Vector{Int64}}([31, 203]) 11038476695863080329 with 2 children and no parent.

julia> using AbstractTrees

julia> print_tree(T, maxdepth=16)
[31, 203]
├─ [29, 171]
│  ├─ [27, 141]
│  │  ├─ [23, 109]
│  │  │  ├─ [21, 85]
│  │  │  │  ├─ [7, 17]
│  │  │  │  │  ├─ [3, 5]
│  │  │  │  │  │  ├─ [1, 1]
│  │  │  │  │  │  └─ [1, 1]
│  │  │  │  │  └─ [3, 5]
│  │  │  │  │     ├─ [1, 1]
│  │  │  │  │     └─ [1, 1]
│  │  │  │  └─ [13, 47]
│  │  │  │     ├─ [11, 33]
│  │  │  │     │  ├─ [5, 11]
│  │  │  │     │  │  ├─ [1, 1]
│  │  │  │     │  │  └─ [3, 5]
│  │  │  │     │  │     ├─ [1, 1]
│  │  │  │     │  │     └─ [1, 1]
│  │  │  │     │  └─ [5, 11]
│  │  │  │     │     ├─ [1, 1]
│  │  │  │     │     └─ [3, 5]
│  │  │  │     │        ├─ [1, 1]
│  │  │  │     │        └─ [1, 1]
│  │  │  │     └─ [1, 1]
│  │  │  └─ [1, 1]
│  │  └─ [3, 5]
│  │     ├─ [1, 1]
│  │     └─ [1, 1]
│  └─ [1, 1]
└─ [1, 1]