Cannot do ML estimation on BD trees

[richelbilderbeek]

The maximum likelihood estimation of the MBD package fails on (Single)BD trees.

# Simulate a BD tree
set.seed(10)
lambda <- 0.3
mu <- 0.1
phylogeny <- mbd_sim_checked(
  mbd_params = create_mbd_params(lambda = lambda, mu = mu, nu = 0.0, q = 0.0),
  crown_age = 2,
  conditioned_on = "non_extinction"
)$tes

# Maximum likelihood of BD tree
ml_est <- mbd_calc_max_lik(
  branching_times = ape::branching.times(phylogeny),
  init_param_values = create_mbd_params(
    lambda = lambda,
    mu = mu,
    nu = 0.0,
    q = 0.0
  ),
  estimated_params = create_mbd_params_selector(lambda = TRUE, mu = TRUE),
  fixed_params = create_mbd_params_selector(nu = TRUE, q = TRUE),
  init_n_species = 2,
  n_missing_species = 0,
  conditioned_on = "non_extinction"
)
expect_true(ml_est$lambda >= 0)

[richelbilderbeek]

Marked test with

skip("Cannot do ML estimation on BD trees, Issue #4, #4")

[Giappo]

Please give a look inside mbd_loglik on the error triggered by condition2. I suppose you get the error because you are using q=0 which is not allowed by the "safety_threshold". In general q==0 could create problems. However if you want to test under BD conditions, any value of q should work, as long as nu==0 (if the environmental change never triggers, the probability of multiple speciation is not influent).

[richelbilderbeek]

Thanks! I'm wondering if some elegant rewrite could allow for q == 0. Combining the code with the mathematical derivation may help me succeed :+1:

[Giappo]

i think you can add an if statement like if (nu == 0) then (whatever with q)

[richelbilderbeek]

Note to self pmb_loglik could handle every combination of q == 0 and nu == 0. This was disabled by an if-statement, but dropping that restriction just worked. For mbd_loglik this restriction cannot be dropped, due to log(pc) being NaN, where pc is created by calc_cond_prob.