Formulas in aphylo — aphylo-model • aphylo

This function the the workhorse behind the likelihood function. It creates arbitrary models by modifying the call to LogLike() function according to what the user specifies as model.

eta(..., env)

psi(..., env)

Pi(..., env)

mu_d(..., env)

mu_s(..., env)

aphylo_formula(fm, params, priors, env = parent.frame())

Arguments

...: Either 0, 1 or both. Depending on the parameter, the index of the model parameter that will be set as fixed.
env: Environment (not to be called by the user).
fm: A formula. Model of the type <aphylo-object> ~ <parameters> (see examples).
params: Numeric vector with model parameters.
priors: (optional) A function. Prior for the model.

Value

A list with the following elements:

fun A function. The log-likelihood function.
fixed Logical vector.

Examples

set.seed(12)
x <- raphylo(10)

# Baseline model
aphylo_formula(x ~ mu_d)
#> $model
#> x ~ mu_d
#> <environment: 0x558979b88d50>
#> 
#> $dat
#> 
#> Phylogenetic tree with 10 tips and 9 internal nodes.
#> 
#> Tip labels:
#>   1, 2, 3, 4, 5, 6, ...
#> Node labels:
#>   11, 12, 13, 14, 15, 16, ...
#> 
#> Rooted; no branch lengths.
#> 
#>  Tip (leafs) annotations:
#>   fun0000
#> 1       1
#> 2       1
#> 3       1
#> 4       1
#> 5       1
#> 6       1
#> 
#> ...(4 obs. omitted)...
#> 
#> 
#>  Internal node annotations:
#>   fun0000
#> 1       1
#> 2       1
#> 3       1
#> 4       1
#> 5       1
#> 6       1
#> 
#> ...(3 obs. omitted)...
#> 
#> 
#> $params
#> mu_d0 mu_d1 
#>   0.9   0.5 
#> 
#> $fixed
#> mu_d0 mu_d1 
#> FALSE FALSE 
#> 
#> $fun
#> function (p, dat, priors, verb_ans = FALSE) 
#> {
#>     ans <- LogLike(tree = dat, psi = c(0, 0), mu_d = p[c("mu_d0", 
#>         "mu_d1")], mu_s = p[c("mu_d0", "mu_d1")], eta = -c(1, 
#>         1), Pi = -1, verb_ans = verb_ans)
#>     ans$ll <- ans$ll + sum(log(priors(p)))
#>     if (!is.finite(ans$ll)) {
#>         ans$ll <- -.Machine$double.xmax * 1e-10
#>     }
#>     if (verb_ans) {
#>         ans
#>     }
#>     else {
#>         ans$ll
#>     }
#> }
#> <bytecode: 0x55897e0ced60>
#> <environment: 0x558979d234d8>
#> 

# Mislabeling probabilities
aphylo_formula(x ~ mu_d + psi)
#> $model
#> x ~ mu_d + psi
#> <environment: 0x558979b88d50>
#> 
#> $dat
#> 
#> Phylogenetic tree with 10 tips and 9 internal nodes.
#> 
#> Tip labels:
#>   1, 2, 3, 4, 5, 6, ...
#> Node labels:
#>   11, 12, 13, 14, 15, 16, ...
#> 
#> Rooted; no branch lengths.
#> 
#>  Tip (leafs) annotations:
#>   fun0000
#> 1       1
#> 2       1
#> 3       1
#> 4       1
#> 5       1
#> 6       1
#> 
#> ...(4 obs. omitted)...
#> 
#> 
#>  Internal node annotations:
#>   fun0000
#> 1       1
#> 2       1
#> 3       1
#> 4       1
#> 5       1
#> 6       1
#> 
#> ...(3 obs. omitted)...
#> 
#> 
#> $params
#>  psi0  psi1 mu_d0 mu_d1 
#>  0.10  0.05  0.90  0.50 
#> 
#> $fixed
#>  psi0  psi1 mu_d0 mu_d1 
#> FALSE FALSE FALSE FALSE 
#> 
#> $fun
#> function (p, dat, priors, verb_ans = FALSE) 
#> {
#>     ans <- LogLike(tree = dat, psi = c(p["psi0"], p["psi1"]), 
#>         mu_d = p[c("mu_d0", "mu_d1")], mu_s = p[c("mu_d0", "mu_d1")], 
#>         eta = -c(1, 1), Pi = -1, verb_ans = verb_ans)
#>     ans$ll <- ans$ll + sum(log(priors(p)))
#>     if (!is.finite(ans$ll)) {
#>         ans$ll <- -.Machine$double.xmax * 1e-10
#>     }
#>     if (verb_ans) {
#>         ans
#>     }
#>     else {
#>         ans$ll
#>     }
#> }
#> <environment: 0x558979dbf2d0>
#> 

# Different probabilities for speciation and duplication node
# (only works if you have both types)
aphylo_formula(x ~ mu_d + mu_s + psi)
#> $model
#> x ~ mu_d + mu_s + psi
#> <environment: 0x558979b88d50>
#> 
#> $dat
#> 
#> Phylogenetic tree with 10 tips and 9 internal nodes.
#> 
#> Tip labels:
#>   1, 2, 3, 4, 5, 6, ...
#> Node labels:
#>   11, 12, 13, 14, 15, 16, ...
#> 
#> Rooted; no branch lengths.
#> 
#>  Tip (leafs) annotations:
#>   fun0000
#> 1       1
#> 2       1
#> 3       1
#> 4       1
#> 5       1
#> 6       1
#> 
#> ...(4 obs. omitted)...
#> 
#> 
#>  Internal node annotations:
#>   fun0000
#> 1       1
#> 2       1
#> 3       1
#> 4       1
#> 5       1
#> 6       1
#> 
#> ...(3 obs. omitted)...
#> 
#> 
#> $params
#>  psi0  psi1 mu_d0 mu_d1 mu_s0 mu_s1 
#>  0.10  0.05  0.90  0.50  0.10  0.05 
#> 
#> $fixed
#>  psi0  psi1 mu_d0 mu_d1 mu_s0 mu_s1 
#> FALSE FALSE FALSE FALSE FALSE FALSE 
#> 
#> $fun
#> function (p, dat, priors, verb_ans = FALSE) 
#> {
#>     ans <- LogLike(tree = dat, psi = c(p["psi0"], p["psi1"]), 
#>         mu_d = p[c("mu_d0", "mu_d1")], mu_s = c(p["mu_s0"], p["mu_s1"]), 
#>         eta = -c(1, 1), Pi = -1, verb_ans = verb_ans)
#>     ans$ll <- ans$ll + sum(log(priors(p)))
#>     if (!is.finite(ans$ll)) {
#>         ans$ll <- -.Machine$double.xmax * 1e-10
#>     }
#>     if (verb_ans) {
#>         ans
#>     }
#>     else {
#>         ans$ll
#>     }
#> }
#> <environment: 0x558979e6f630>
#> 

# Mislabeling probabilities and etas(fixed)
aphylo_formula(x ~ mu_d + psi + eta(0, 1))
#> $model
#> x ~ mu_d + psi + eta(0, 1)
#> <environment: 0x558979b88d50>
#> 
#> $dat
#> 
#> Phylogenetic tree with 10 tips and 9 internal nodes.
#> 
#> Tip labels:
#>   1, 2, 3, 4, 5, 6, ...
#> Node labels:
#>   11, 12, 13, 14, 15, 16, ...
#> 
#> Rooted; no branch lengths.
#> 
#>  Tip (leafs) annotations:
#>   fun0000
#> 1       1
#> 2       1
#> 3       1
#> 4       1
#> 5       1
#> 6       1
#> 
#> ...(4 obs. omitted)...
#> 
#> 
#>  Internal node annotations:
#>   fun0000
#> 1       1
#> 2       1
#> 3       1
#> 4       1
#> 5       1
#> 6       1
#> 
#> ...(3 obs. omitted)...
#> 
#> 
#> $params
#>  psi0  psi1 mu_d0 mu_d1  eta0  eta1 
#>  0.10  0.05  0.90  0.50  1.00  1.00 
#> 
#> $fixed
#>  psi0  psi1 mu_d0 mu_d1  eta0  eta1 
#> FALSE FALSE FALSE FALSE  TRUE  TRUE 
#> 
#> $fun
#> function (p, dat, priors, verb_ans = FALSE) 
#> {
#>     ans <- LogLike(tree = dat, psi = c(p["psi0"], p["psi1"]), 
#>         mu_d = p[c("mu_d0", "mu_d1")], mu_s = p[c("mu_d0", "mu_d1")], 
#>         eta = c(p["eta0"], p["eta1"]), Pi = -1, verb_ans = verb_ans)
#>     ans$ll <- ans$ll + sum(log(priors(p)))
#>     if (!is.finite(ans$ll)) {
#>         ans$ll <- -.Machine$double.xmax * 1e-10
#>     }
#>     if (verb_ans) {
#>         ans
#>     }
#>     else {
#>         ans$ll
#>     }
#> }
#> <environment: 0x558979f1a0e8>
#> 

# Mislabeling probabilities and Pi 
aphylo_formula(x ~ mu_d + psi + Pi)
#> $model
#> x ~ mu_d + psi + Pi
#> <environment: 0x558979b88d50>
#> 
#> $dat
#> 
#> Phylogenetic tree with 10 tips and 9 internal nodes.
#> 
#> Tip labels:
#>   1, 2, 3, 4, 5, 6, ...
#> Node labels:
#>   11, 12, 13, 14, 15, 16, ...
#> 
#> Rooted; no branch lengths.
#> 
#>  Tip (leafs) annotations:
#>   fun0000
#> 1       1
#> 2       1
#> 3       1
#> 4       1
#> 5       1
#> 6       1
#> 
#> ...(4 obs. omitted)...
#> 
#> 
#>  Internal node annotations:
#>   fun0000
#> 1       1
#> 2       1
#> 3       1
#> 4       1
#> 5       1
#> 6       1
#> 
#> ...(3 obs. omitted)...
#> 
#> 
#> $params
#>  psi0  psi1 mu_d0 mu_d1    Pi 
#>  0.10  0.05  0.90  0.50  0.50 
#> 
#> $fixed
#>  psi0  psi1 mu_d0 mu_d1    Pi 
#> FALSE FALSE FALSE FALSE FALSE 
#> 
#> $fun
#> function (p, dat, priors, verb_ans = FALSE) 
#> {
#>     ans <- LogLike(tree = dat, psi = c(p["psi0"], p["psi1"]), 
#>         mu_d = p[c("mu_d0", "mu_d1")], mu_s = p[c("mu_d0", "mu_d1")], 
#>         eta = -c(1, 1), Pi = p["Pi"], verb_ans = verb_ans)
#>     ans$ll <- ans$ll + sum(log(priors(p)))
#>     if (!is.finite(ans$ll)) {
#>         ans$ll <- -.Machine$double.xmax * 1e-10
#>     }
#>     if (verb_ans) {
#>         ans
#>     }
#>     else {
#>         ans$ll
#>     }
#> }
#> <environment: 0x558979fc6660>
#>