Mainland extinction correct?

[richelbilderbeek]

Dear DAISIEmainland maintainer,

Thanks for DAISIEmainland and its well-written and well-tested code! I've already mentioned that I missed a function to plot mainland and island, so -up till now- I wrote a function to plot the mainland (it is on branch richel). Within these plots, I do see things I do not understand. Please correct me.

So, here is my first test, in test-plot_mainland.R:

test_that("minimal example use", {
  set.seed(
    1,
    kind = "Mersenne-Twister",
    normal.kind = "Inversion",
    sample.kind = "Rejection"
  )
  mainland <- DAISIEmainland:::sim_mainland(
    total_time = 1,
    m = 2,
    mainland_ex = 1)

  plot_mainland(mainland)
})

This is the code and plot:

This is great! The number of species in all clades in nicely preserved (I first thought this was a mistake, but this 'duo pattern' always happens when there is only 1 clade).

Now to a more interesting example, with only a different seed and extinction rate:

test_that("interesting example", {
  set.seed(
    4,
    kind = "Mersenne-Twister",
    normal.kind = "Inversion",
    sample.kind = "Rejection"
  )
  mainland <- DAISIEmainland:::sim_mainland(
    total_time = 1,
    m = 4,
    mainland_ex = 1)

  plot_mainland(mainland)
})

Now, this is something I'd like to be checked upon:

At t = 0.04290015 (see mainland printout below), in clades 3 and 4, the same extinction seems to happen and this results in a nice duo in clade 4. Should it be possible that two species from different clades go extinct at exactly the same time?

The answer: yes! When species 15 (in clade 3) goes extinct, it lets the ancestral species in clade 4 speciate!

> mainland
[[1]]
  spec_id main_anc_id spec_type branch_code  branch_t spec_origin_t spec_ex_t
1       1           1         E           A       NaN     0.0000000 0.3282561
2      11           1         E          AA 0.3282561     0.3282561 0.4520710
3      12           1         E          AB 0.3282561     0.3282561 0.4520710
4      13           1         E         ABA 0.4520710     0.4520710 0.7088819
5      14           1         E         ABB 0.4520710     0.4520710 0.8010316
6      15           1         E        ABAA 0.7088819     0.7088819 0.7237706
7      16           1         C        ABAB 0.7088819     0.7088819 1.0000000

[[2]]
  spec_id main_anc_id spec_type branch_code  branch_t spec_origin_t spec_ex_t
1       2           2         E           A       NaN     0.0000000 0.2562369
2       9           2         E          AA 0.2562369     0.2562369 0.7237706
3      10           2         E          AB 0.2562369     0.2562369 0.3282561
4      17           2         E         AAA 0.7237706     0.7237706 0.8010316
5      18           2         C         AAB 0.7237706     0.7237706 1.0000000
6      19           2         C        AAAA 0.8010316     0.8010316 1.0000000
7      20           2         C        AAAB 0.8010316     0.8010316 1.0000000

[[3]]
  spec_id main_anc_id spec_type branch_code branch_t spec_origin_t  spec_ex_t
1       3           3         E           A      NaN             0 0.04290015

[[4]]
  spec_id main_anc_id spec_type branch_code   branch_t spec_origin_t  spec_ex_t
1       4           4         E           A        NaN    0.00000000 0.04290015
2       5           4         E          AA 0.04290015    0.04290015 0.24356191
3       6           4         E          AB 0.04290015    0.04290015 0.24356191
4       7           4         E         ABA 0.24356191    0.24356191 0.25623688
5       8           4         E         ABB 0.24356191    0.24356191 0.70888185```

[richelbilderbeek]

&TLDR: everything works as expected :-)