Using secsse with complete phylogenies (with extinction)

Pedro Santos Neves

2025-07-24

Introduction

Most current studies of evolutionary dynamics make use of molecular phylogenies, which, for most groups, contain only information on extant species. However, when data on extinct species is available, usually through the presence of fossil data, we can use complete trees. Thus, we can leverage the data from extinct lineages for maximum-likelihood estimation with secsse.

Note that here “complete tree” should not be taken as a complete sampling fraction, that is, all known species being present in the phylogeny and there being no missing data, but rather the assumption that all currently extinct species are included. This follows the nomenclature of Nee et al. (1994), who also coined the term “reconstructed tree” for phylogenies for which there is no information on extinct lineages.

Set-up

Like all ML analyses with secsse, we first need a few things to start with, starting with a dated phylogeny. For the purpose of this vignette, we are going to simulate phylogenies with secsse_sim(). We will simulate a reconstructed and a complete version of the same tree under the CR model.

In order to simulate the trees, we need to specify the model and set starting parameters. Here we simulate a similar constant rate (CR) example to that of the Simulating with secsse vignette. For more details on this model and full details of the functionality of secsse_sim(), see vignette("sim_with_secsse", package = "secsse").

library(secsse)

spec_matrix <- c()
spec_matrix <- rbind(spec_matrix, c(0, 0, 0, 1))
spec_matrix <- rbind(spec_matrix, c(1, 1, 1, 1))
lambda_list <- secsse::create_lambda_list(state_names = c(0, 1),
                                          num_concealed_states = 2,
                                          transition_matrix = spec_matrix,
                                          model = "CR")

mu_vector <- secsse::create_mu_vector(state_names = c(0, 1),
                                      num_concealed_states = 2,
                                      model = "CR",
                                      lambda_list = lambda_list)

shift_matrix <- c()
shift_matrix <- rbind(shift_matrix, c(0, 1, 3))
shift_matrix <- rbind(shift_matrix, c(1, 0, 4))

q_matrix <- secsse::create_q_matrix(state_names = c(0, 1),
                                    num_concealed_states = 2,
                                    shift_matrix = shift_matrix,
                                    diff.conceal = FALSE)

# Set-up starting parameters
speciation_rate <- 0.8
extinction_rate <- 0.2
q_01 <- 0.1
q_10 <- 0.1
used_params <- c(speciation_rate, extinction_rate, q_01, q_10)

sim_lambda_list <- secsse::fill_in(lambda_list, used_params)
sim_mu_vector   <- secsse::fill_in(mu_vector, used_params)
sim_q_matrix    <- secsse::fill_in(q_matrix, used_params)

# Simulate and plot the tree

sim_tree_complete <- secsse::secsse_sim(lambdas = sim_lambda_list,
                                        mus = sim_mu_vector,
                                        qs = sim_q_matrix,
                                        crown_age = 5,
                                        num_concealed_states = 2,
                                        seed = 40,
                                        drop_extinct = FALSE)

if (requireNamespace("diversitree")) {
  traits_for_plot_complete <- data.frame(
    trait = as.numeric(sim_tree_complete$obs_traits),
    row.names = sim_tree_complete$phy$tip.label
  )
  diversitree::trait.plot(tree = sim_tree_complete$phy,
                          dat = traits_for_plot_complete,
                          cols = list("trait" = c("blue", "red")),
                          type = "p")
} else {
  plot(sim_tree_complete$phy)
}
#> Loading required namespace: diversitree

Fitting the model

Finally, we run secsse_ml() on our complete tree, much in the same way as we would for one with extant species. However, this time we make sure to set the is_complete_tree argument to TRUE (defaults to FALSE if omitted). This enables secsse to use the information present in extinct lineages.

idparsopt <- 1:4 # our maximum rate parameter was 4 -> We are keeping
# concealed and examined traits the same for the MLE.
idparsfix <- c(0) # we want to keep all zeros at zero
initparsopt <- rep(0.1, 4)
initparsfix <- c(0.0) # all zeros remain at zero.
sampling_fraction <- c(1, 1)

idparslist <- list()
idparslist[[1]] <- lambda_list
idparslist[[2]] <- mu_vector
idparslist[[3]] <- q_matrix

complete_tree_ml_CR <- secsse_ml(phy = sim_tree_complete$phy,
                                 traits = sim_tree_complete$obs_traits,
                                 num_concealed_states = 2,
                                 idparslist = idparslist,
                                 idparsopt = idparsopt,
                                 initparsopt = initparsopt,
                                 idparsfix = idparsfix,
                                 parsfix = initparsfix,
                                 sampling_fraction = sampling_fraction,
                                 verbose = FALSE)
#> Note: you set some transitions as impossible to happen.
#> 1 0.1 0.1 0.1 0.1 -333.41949072334 initial 
#> 2 0.105 0.1 0.1 0.1 -328.399396278715 reflect 
#> 3 0.105 0.1 0.1 0.1 -328.399396278715 reflect 
#> 4 0.10846394984326 0.0944288126055149 0.0990675332014694 0.10846394984326 -324.265242339137 expand 
#> 5 0.114897886062343 0.0902873160476524 0.0983692199082245 0.101011959521619 -318.155872449178 expand 
#> 6 0.114897886062343 0.0902873160476524 0.0983692199082245 0.101011959521619 -318.155872449178 reflect 
#> 7 0.124172002998582 0.0905356861233287 0.101313961616394 0.102660370839572 -310.361759023057 expand 
#> 8 0.136331330154724 0.0787355101510194 0.0948548055460873 0.112245294610444 -299.669043268782 expand 
#> 9 0.136331330154724 0.0787355101510194 0.0948548055460873 0.112245294610444 -299.669043268782 reflect 
#> 10 0.157463962842126 0.0655976490744996 0.103698336063287 0.104834174336294 -284.755391015959 expand 
#> 11 0.188580384493981 0.0584948324319094 0.0996649403718381 0.11398400869047 -267.200567892714 expand 
#> 12 0.188580384493981 0.0584948324319094 0.0996649403718381 0.11398400869047 -267.200567892714 reflect 
#> 13 0.232811307810817 0.0269459250512547 0.104513689147747 0.130943534064202 -246.325665194653 expand 
#> 14 0.250333207168013 0.0242967293163897 0.106866374342769 0.119320701155624 -240.34682488793 reflect 
#> 15 0.275886868150326 0.0164978330464061 0.0995375741889073 0.134415859843176 -232.074212604131 reflect 
#> 16 0.290567786992314 0.0104768915751455 0.109885889970718 0.13578248786944 -227.902107778429 reflect 
#> 17 0.290567786992314 0.0104768915751455 0.109885889970718 0.13578248786944 -227.902107778429 contract inside 
#> 18 0.302225873216359 0.00107913758124338 0.107970254509982 0.138301143863997 -224.484253738002 reflect 
#> 19 0.302225873216359 0.00107913758124338 0.107970254509982 0.138301143863997 -224.484253738002 contract inside 
#> 20 0.302225873216359 0.00107913758124338 0.107970254509982 0.138301143863997 -224.484253738002 contract inside 
#> 21 0.312285510183205 0.00318560730171753 0.106518346468716 0.136352050445652 -222.2740881016 reflect 
#> 22 0.312285510183205 0.00318560730171753 0.106518346468716 0.136352050445652 -222.2740881016 contract inside 
#> 23 0.312285510183205 0.00318560730171753 0.106518346468716 0.136352050445652 -222.2740881016 contract inside 
#> 24 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -221.358334181637 reflect 
#> 25 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -221.358334181637 contract inside 
#> 26 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -221.358334181637 contract inside 
#> 27 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -221.358334181637 contract inside 
#> 28 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -221.358334181637 contract inside 
#> 29 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -221.358334181637 contract inside 
#> 30 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -221.358334181637 reflect 
#> 31 0.321279685841566 0.00148420430175912 0.106255055896489 0.139103791496767 -220.208484370835 reflect 
#> 32 0.327951378167453 0.00116396241305509 0.106520668879064 0.139522975244395 -218.792001269951 reflect 
#> 33 0.328159974506911 0.000744716688942316 0.107153785306222 0.142586871134893 -218.736375957232 reflect 
#> 34 0.328159974506911 0.000744716688942316 0.107153785306222 0.142586871134893 -218.736375957232 contract inside 
#> 35 0.341455582607457 0.00404681386379073 0.104381221672472 0.138214985108092 -216.224415091154 expand 
#> 36 0.341455582607457 0.00404681386379073 0.104381221672472 0.138214985108092 -216.224415091154 reflect 
#> 37 0.361306082996556 0.00210566130103047 0.105020910393752 0.142307910206395 -212.498054458597 expand 
#> 38 0.361306082996556 0.00210566130103047 0.105020910393752 0.142307910206395 -212.498054458597 reflect 
#> 39 0.395988449433568 0.0057685064665477 0.100079074063907 0.137159385636665 -207.061282116041 expand 
#> 40 0.395988449433568 0.0057685064665477 0.100079074063907 0.137159385636665 -207.061282116041 reflect 
#> 41 0.449547208116093 0.00469530325609993 0.0995992258414071 0.144521055959469 -200.134693372523 expand 
#> 42 0.449547208116093 0.00469530325609993 0.0995992258414071 0.144521055959469 -200.134693372523 reflect 
#> 43 0.547809355014861 0.013735148623196 0.0901324106268667 0.135999345118161 -191.634486188102 expand 
#> 44 0.547809355014861 0.013735148623196 0.0901324106268667 0.135999345118161 -191.634486188102 reflect 
#> 45 0.725947805699673 0.0141854163815616 0.0857802257181384 0.145022748857269 -184.7165341319 expand 
#> 46 0.725947805699673 0.0141854163815616 0.0857802257181384 0.145022748857269 -184.7165341319 reflect 
#> 47 0.836859697016871 0.0208632470773777 0.0775784053345886 0.13686049855624 -184.130865310278 reflect 
#> 48 0.836859697016871 0.0208632470773777 0.0775784053345886 0.13686049855624 -184.130865310278 reflect 
#> 49 0.836859697016871 0.0208632470773777 0.0775784053345886 0.13686049855624 -184.130865310278 contract outside 
#> 50 0.836859697016871 0.0208632470773777 0.0775784053345886 0.13686049855624 -184.130865310278 contract inside 
#> 51 0.836859697016871 0.0208632470773777 0.0775784053345886 0.13686049855624 -184.130865310278 contract inside 
#> 52 0.836859697016871 0.0208632470773777 0.0775784053345886 0.13686049855624 -184.130865310278 contract inside 
#> 53 0.836859697016871 0.0208632470773777 0.0775784053345886 0.13686049855624 -184.130865310278 contract inside 
#> 54 0.791695761391743 0.0178010678566669 0.0811135164792597 0.141550153685247 -184.10592505601 contract inside 
#> 55 0.788978648668157 0.018793490591944 0.0805542323109085 0.139631107972986 -184.100130479146 contract outside 
#> 56 0.793687597982675 0.0179916159930553 0.0811283837425861 0.141390063708646 -184.099479144292 contract inside 
#> 57 0.822723340912633 0.0194384590416029 0.0795761713600517 0.140660439989321 -184.097339437644 contract inside 
#> 58 0.817818871624196 0.019683129812061 0.0790834512702414 0.138830548426376 -184.07701462465 contract inside 
#> 59 0.817818871624196 0.019683129812061 0.0790834512702414 0.138830548426376 -184.07701462465 reflect 
#> 60 0.801133909653241 0.0190547783631467 0.0801324190153936 0.139763557018266 -184.07481776991 contract inside 
#> 61 0.801133909653241 0.0190547783631467 0.0801324190153936 0.139763557018266 -184.07481776991 contract inside 
#> 62 0.799007648687733 0.0194002244963925 0.0797077257740846 0.138211891635282 -184.070216853191 reflect 
#> 63 0.81439899381252 0.0203606167427226 0.0786970267714734 0.137321353393389 -184.061854420671 reflect 
#> 64 0.81439899381252 0.0203606167427226 0.0786970267714734 0.137321353393389 -184.061854420671 reflect 
#> 65 0.81439899381252 0.0203606167427226 0.0786970267714734 0.137321353393389 -184.061854420671 contract inside 
#> 66 0.81439899381252 0.0203606167427226 0.0786970267714734 0.137321353393389 -184.061854420671 reflect 
#> 67 0.81439899381252 0.0203606167427226 0.0786970267714734 0.137321353393389 -184.061854420671 reflect 
#> 68 0.81439899381252 0.0203606167427226 0.0786970267714734 0.137321353393389 -184.061854420671 contract inside 
#> 69 0.805565812585464 0.02009247263236 0.0783726703278158 0.13571930509537 -184.056056935955 expand 
#> 70 0.813717263007625 0.0215241196832835 0.0773641149124775 0.133516567856987 -184.048162763931 expand 
#> 71 0.813717263007625 0.0215241196832835 0.0773641149124775 0.133516567856987 -184.048162763931 reflect 
#> 72 0.812622242256909 0.0220682745808221 0.0761445548159273 0.130683566680483 -184.045769766894 expand 
#> 73 0.812622242256909 0.0220682745808221 0.0761445548159273 0.130683566680483 -184.045769766894 reflect 
#> 74 0.812622242256909 0.0220682745808221 0.0761445548159273 0.130683566680483 -184.045769766894 reflect 
#> 75 0.812622242256909 0.0220682745808221 0.0761445548159273 0.130683566680483 -184.045769766894 contract outside 
#> 76 0.80812764845636 0.0218968818918979 0.0763335132221221 0.130549507309141 -184.044891401942 contract inside 
#> 77 0.808577842499566 0.0217245482811683 0.0769830090762486 0.132189765232994 -184.043767064511 reflect 
#> 78 0.808577842499566 0.0217245482811683 0.0769830090762486 0.132189765232994 -184.043767064511 contract outside 
#> 79 0.808577842499566 0.0217245482811683 0.0769830090762486 0.132189765232994 -184.043767064511 contract inside 
#> 80 0.808577842499566 0.0217245482811683 0.0769830090762486 0.132189765232994 -184.043767064511 reflect 
#> 81 0.808577842499566 0.0217245482811683 0.0769830090762486 0.132189765232994 -184.043767064511 reflect 
#> 82 0.809692102519732 0.0221932031604373 0.0769745696378608 0.131972741449657 -184.040755185026 expand 
#> 83 0.809692102519732 0.0221932031604373 0.0769745696378608 0.131972741449657 -184.040755185026 contract inside 
#> 84 0.809692102519732 0.0221932031604373 0.0769745696378608 0.131972741449657 -184.040755185026 reflect 
#> 85 0.809692102519732 0.0221932031604373 0.0769745696378608 0.131972741449657 -184.040755185026 reflect 
#> 86 0.813809835121397 0.0232544361830925 0.0759273082591948 0.129617624795755 -184.039847298242 expand 
#> 87 0.813996173264448 0.0231516449535424 0.0768921445660287 0.131954349565469 -184.036408357894 expand 
#> 88 0.813996173264448 0.0231516449535424 0.0768921445660287 0.131954349565469 -184.036408357894 reflect 
#> 89 0.817914294027173 0.0241060274927162 0.0769219683520157 0.131734336734639 -184.034266477645 expand 
#> 90 0.817914294027173 0.0241060274927162 0.0769219683520157 0.131734336734639 -184.034266477645 reflect 
#> 91 0.818449776170596 0.0250771152115688 0.0774522345543237 0.132645692091794 -184.026143188344 expand 
#> 92 0.818449776170596 0.0250771152115688 0.0774522345543237 0.132645692091794 -184.026143188344 reflect 
#> 93 0.818449776170596 0.0250771152115688 0.0774522345543237 0.132645692091794 -184.026143188344 reflect 
#> 94 0.818449776170596 0.0250771152115688 0.0774522345543237 0.132645692091794 -184.026143188344 reflect 
#> 95 0.821410264873652 0.0259187799894594 0.0788123189565693 0.135206847991287 -184.024855716769 expand 
#> 96 0.82548464435546 0.028626783808919 0.0781931289861355 0.133102049113227 -184.010807022587 expand 
#> 97 0.82548464435546 0.028626783808919 0.0781931289861355 0.133102049113227 -184.010807022587 reflect 
#> 98 0.82548464435546 0.028626783808919 0.0781931289861355 0.133102049113227 -184.010807022587 reflect 
#> 99 0.827936648906467 0.0310196164342713 0.0817411565908304 0.1395244809963 -184.006419092006 expand 
#> 100 0.826876397378651 0.0331672228076757 0.0812056026024146 0.137605238803379 -183.979922184496 expand 
#> 101 0.826876397378651 0.0331672228076757 0.0812056026024146 0.137605238803379 -183.979922184496 reflect 
#> 102 0.826876397378651 0.0331672228076757 0.0812056026024146 0.137605238803379 -183.979922184496 reflect 
#> 103 0.835624062578356 0.0419655963912653 0.0869624349196876 0.14608127948256 -183.970021089187 expand 
#> 104 0.836462266334855 0.0453741622443226 0.0837187047688226 0.138017559493735 -183.905061680971 expand 
#> 105 0.836462266334855 0.0453741622443226 0.0837187047688226 0.138017559493735 -183.905061680971 reflect 
#> 106 0.836462266334855 0.0453741622443226 0.0837187047688226 0.138017559493735 -183.905061680971 reflect 
#> 107 0.847486540015565 0.0659863345816972 0.0946625269733423 0.152274207271603 -183.866869762432 expand 
#> 108 0.838676067420241 0.06584295898488 0.0887708438868339 0.139816699462145 -183.759895461832 expand 
#> 109 0.838676067420241 0.06584295898488 0.0887708438868339 0.139816699462145 -183.759895461832 reflect 
#> 110 0.857926860689731 0.0918215096957005 0.0942573085624215 0.141176879149378 -183.626146310481 expand 
#> 111 0.857926860689731 0.0918215096957005 0.0942573085624215 0.141176879149378 -183.626146310481 reflect 
#> 112 0.866919614276863 0.121173316293174 0.0952879149046811 0.133166648287249 -183.436104233893 expand 
#> 113 0.866919614276863 0.121173316293174 0.0952879149046811 0.133166648287249 -183.436104233893 reflect 
#> 114 0.908971127793232 0.206280732222862 0.113832418015663 0.141803726909618 -183.122354997902 expand 
#> 115 0.908971127793232 0.206280732222862 0.113832418015663 0.141803726909618 -183.122354997902 reflect 
#> 116 0.921543076750175 0.29078278713034 0.113439366790698 0.116842085015816 -183.012002848422 expand 
#> 117 0.972825552572185 0.372598561082764 0.122859072769269 0.116598438429242 -182.888254101464 expand 
#> 118 0.972825552572185 0.372598561082764 0.122859072769269 0.116598438429242 -182.888254101464 reflect 
#> 119 0.972825552572185 0.372598561082764 0.122859072769269 0.116598438429242 -182.888254101464 contract outside 
#> 120 0.972825552572185 0.372598561082764 0.122859072769269 0.116598438429242 -182.888254101464 contract inside 
#> 121 0.972825552572185 0.372598561082764 0.122859072769269 0.116598438429242 -182.888254101464 reflect 
#> 122 1.01829134998216 0.416421931783558 0.139015136322224 0.138785131933327 -182.774732109189 expand 
#> 123 1.01829134998216 0.416421931783558 0.139015136322224 0.138785131933327 -182.774732109189 contract inside 
#> 124 1.01829134998216 0.416421931783558 0.139015136322224 0.138785131933327 -182.774732109189 reflect 
#> 125 1.01642557565997 0.412499567933597 0.132431248073012 0.12806200962451 -182.752320413384 reflect 
#> 126 1.01642557565997 0.412499567933597 0.132431248073012 0.12806200962451 -182.752320413384 contract inside 
#> 127 1.01642557565997 0.412499567933597 0.132431248073012 0.12806200962451 -182.752320413384 reflect 
#> 128 1.01642557565997 0.412499567933597 0.132431248073012 0.12806200962451 -182.752320413384 reflect 
#> 129 1.01642557565997 0.412499567933597 0.132431248073012 0.12806200962451 -182.752320413384 reflect 
#> 130 1.01642557565997 0.412499567933597 0.132431248073012 0.12806200962451 -182.752320413384 contract inside 
#> 131 1.01642557565997 0.412499567933597 0.132431248073012 0.12806200962451 -182.752320413384 contract inside 
#> 132 1.01642557565997 0.412499567933597 0.132431248073012 0.12806200962451 -182.752320413384 contract inside 
#> 133 1.01642557565997 0.412499567933597 0.132431248073012 0.12806200962451 -182.752320413384 reflect 
#> 134 1.01642557565997 0.412499567933597 0.132431248073012 0.12806200962451 -182.752320413384 reflect 
#> 135 1.02971360979734 0.42389014841986 0.13593259084302 0.132555232597692 -182.746495501793 contract inside 
#> 136 1.02971360979734 0.42389014841986 0.13593259084302 0.132555232597692 -182.746495501793 contract inside 
#> 137 1.04259686465577 0.432612741133675 0.133643526015259 0.128036147593875 -182.741545685014 expand 
#> 138 1.04259686465577 0.432612741133675 0.133643526015259 0.128036147593875 -182.741545685014 contract inside 
#> 139 1.04259686465577 0.432612741133675 0.133643526015259 0.128036147593875 -182.741545685014 reflect 
#> 140 1.04259686465577 0.432612741133675 0.133643526015259 0.128036147593875 -182.741545685014 contract outside 
#> 141 1.04259686465577 0.432612741133675 0.133643526015259 0.128036147593875 -182.741545685014 contract inside 
#> 142 1.04259686465577 0.432612741133675 0.133643526015259 0.128036147593875 -182.741545685014 contract inside 
#> 143 1.03454274667806 0.430739675406173 0.134791808832817 0.129730069277528 -182.740822793631 reflect 
#> 144 1.03454274667806 0.430739675406173 0.134791808832817 0.129730069277528 -182.740822793631 contract inside 
#> 145 1.04367277437857 0.438715201352587 0.134182160485362 0.127736292529269 -182.740740615173 reflect 
#> 146 1.04136836566154 0.431794035860322 0.133793272118611 0.128566663862566 -182.739279958651 reflect 
#> 147 1.04136836566154 0.431794035860322 0.133793272118611 0.128566663862566 -182.739279958651 reflect 
#> 148 1.03571114200258 0.430300529464868 0.133346626195677 0.127725681425213 -182.73896977368 reflect 
#> 149 1.03571114200258 0.430300529464868 0.133346626195677 0.127725681425213 -182.73896977368 reflect 
#> 150 1.0337091175291 0.426188450991544 0.134030156462182 0.1299001016884 -182.736311372846 expand 
#> 151 1.0337091175291 0.426188450991544 0.134030156462182 0.1299001016884 -182.736311372846 reflect 
#> 152 1.0337091175291 0.426188450991544 0.134030156462182 0.1299001016884 -182.736311372846 reflect 
#> 153 1.0337091175291 0.426188450991544 0.134030156462182 0.1299001016884 -182.736311372846 contract inside 
#> 154 1.0337091175291 0.426188450991544 0.134030156462182 0.1299001016884 -182.736311372846 reflect 
#> 155 1.02788339201769 0.414865444944747 0.132651643887932 0.129905201772973 -182.735664929651 reflect 
#> 156 1.02788339201769 0.414865444944747 0.132651643887932 0.129905201772973 -182.735664929651 contract inside 
#> 157 1.03191721674116 0.424771975015852 0.133544631364935 0.129650479684515 -182.733313020447 expand 
#> 158 1.03191721674116 0.424771975015852 0.133544631364935 0.129650479684515 -182.733313020447 reflect 
#> 159 1.03191721674116 0.424771975015852 0.133544631364935 0.129650479684515 -182.733313020447 reflect 
#> 160 1.0287702646285 0.413389993023254 0.132737925374329 0.131190769053712 -182.731395848171 expand 
#> 161 1.04113889865012 0.435030279483414 0.135563766793147 0.132146625096944 -182.730255235933 expand 
#> 162 1.02688853788288 0.420942659455222 0.134212024488607 0.132519762825846 -182.727760044261 expand 
#> 163 1.03290167738434 0.425268239245973 0.133508321075045 0.131356513478487 -182.719960894202 expand 
#> 164 1.03290167738434 0.425268239245973 0.133508321075045 0.131356513478487 -182.719960894202 reflect 
#> 165 1.03290167738434 0.425268239245973 0.133508321075045 0.131356513478487 -182.719960894202 reflect 
#> 166 1.01603308854146 0.410601863590341 0.132619522902578 0.13444022334681 -182.714520441742 expand 
#> 167 1.0361329235941 0.430678504319105 0.134121501952445 0.135134380301786 -182.699110548723 expand 
#> 168 1.0361329235941 0.430678504319105 0.134121501952445 0.135134380301786 -182.699110548723 reflect 
#> 169 1.00918002214349 0.395860392107273 0.128237221452958 0.133135046855068 -182.690103216137 expand 
#> 170 1.00183495870088 0.399686690503877 0.129507529917724 0.139405608352911 -182.688369978317 expand 
#> 171 1.02450280209356 0.420450546646261 0.128946682008103 0.136912927539471 -182.641373587889 expand 
#> 172 1.02450280209356 0.420450546646261 0.128946682008103 0.136912927539471 -182.641373587889 reflect 
#> 173 1.02450280209356 0.420450546646261 0.128946682008103 0.136912927539471 -182.641373587889 reflect 
#> 174 1.02450280209356 0.420450546646261 0.128946682008103 0.136912927539471 -182.641373587889 reflect 
#> 175 1.00929611989185 0.390410518381561 0.118964915603254 0.14033882851506 -182.575897937221 expand 
#> 176 1.00929611989185 0.390410518381561 0.118964915603254 0.14033882851506 -182.575897937221 reflect 
#> 177 1.03785001139668 0.440483070846318 0.122362217350975 0.143959456193926 -182.532836818509 expand 
#> 178 1.03785001139668 0.440483070846318 0.122362217350975 0.143959456193926 -182.532836818509 reflect 
#> 179 1.05315378389886 0.432600228963474 0.111283113058808 0.14806423891212 -182.47909953758 expand 
#> 180 1.04438857156113 0.42365891695664 0.102229483377311 0.144963169695225 -182.405822685632 expand 
#> 181 1.04438857156113 0.42365891695664 0.102229483377311 0.144963169695225 -182.405822685632 reflect 
#> 182 1.04438857156113 0.42365891695664 0.102229483377311 0.144963169695225 -182.405822685632 reflect 
#> 183 1.04438857156113 0.42365891695664 0.102229483377311 0.144963169695225 -182.405822685632 reflect 
#> 184 1.08224913420321 0.469043207378432 0.0902875132279503 0.151853223245422 -182.368813694856 reflect 
#> 185 1.08224913420321 0.469043207378432 0.0902875132279503 0.151853223245422 -182.368813694856 reflect 
#> 186 1.08224913420321 0.469043207378432 0.0902875132279503 0.151853223245422 -182.368813694856 reflect 
#> 187 1.04351640266287 0.432444217053479 0.0963886146942909 0.136449369447079 -182.331094363225 expand 
#> 188 1.04351640266287 0.432444217053479 0.0963886146942909 0.136449369447079 -182.331094363225 reflect 
#> 189 1.04351640266287 0.432444217053479 0.0963886146942909 0.136449369447079 -182.331094363225 reflect 
#> 190 1.08509589870262 0.460053226640003 0.0795806576431622 0.13974420784199 -182.326138528352 reflect 
#> 191 1.1055195682824 0.532074145400752 0.0989760961742834 0.13885950775035 -182.270142836303 expand 
#> 192 1.0634650549951 0.462601863031994 0.0853954227504189 0.1301074728005 -182.259009659886 reflect 
#> 193 1.0634650549951 0.462601863031994 0.0853954227504189 0.1301074728005 -182.259009659886 contract inside 
#> 194 1.11764804892406 0.524842939930668 0.0773319370201016 0.139056420673453 -182.248622861094 reflect 
#> 195 1.0907633919066 0.527987705434282 0.0929134736999873 0.135423951671827 -182.231021524892 reflect 
#> 196 1.0907633919066 0.527987705434282 0.0929134736999873 0.135423951671827 -182.231021524892 reflect 
#> 197 1.09138811835007 0.508978452403214 0.0759498342365384 0.128159791243063 -182.215777381226 reflect 
#> 198 1.09138811835007 0.508978452403214 0.0759498342365384 0.128159791243063 -182.215777381226 reflect 
#> 199 1.09138811835007 0.508978452403214 0.0759498342365384 0.128159791243063 -182.215777381226 reflect 
#> 200 1.09832526276809 0.554043960053072 0.0850687410039704 0.130624172483883 -182.214940657914 reflect 
#> 201 1.10423062568848 0.53176112497475 0.0780571912559945 0.13618764068599 -182.198198074478 contract inside 
#> 202 1.10423062568848 0.53176112497475 0.0780571912559945 0.13618764068599 -182.198198074478 reflect 
#> 203 1.12754488327672 0.574976934060789 0.08075396092941 0.133560483228174 -182.19816438664 contract inside 
#> 204 1.12754488327672 0.574976934060789 0.08075396092941 0.133560483228174 -182.19816438664 contract inside 
#> 205 1.11948483212239 0.560646982897376 0.0746447498250853 0.131339417181906 -182.19646802865 contract inside 
#> 206 1.10591692867801 0.551988787423421 0.0813366295520978 0.131740909902765 -182.195581132683 contract inside 
#> 207 1.10909842359642 0.543325646791719 0.078252284569194 0.133965190872232 -182.193037888394 contract inside 
#> 208 1.10909842359642 0.543325646791719 0.078252284569194 0.133965190872232 -182.193037888394 contract inside 
#> 209 1.10909842359642 0.543325646791719 0.078252284569194 0.133965190872232 -182.193037888394 contract outside 
#> 210 1.10909842359642 0.543325646791719 0.078252284569194 0.133965190872232 -182.193037888394 contract inside 
#> 211 1.10723682339792 0.548422211551274 0.0793427391477492 0.131951291928255 -182.192416199103 contract inside 
#> 212 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 expand 
#> 213 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 reflect 
#> 214 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 reflect 
#> 215 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 reflect 
#> 216 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 contract inside 
#> 217 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 contract inside 
#> 218 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 reflect 
#> 219 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 reflect 
#> 220 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 reflect 
#> 221 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 contract inside 
#> 222 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 contract inside 
#> 223 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 reflect 
#> 224 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 contract outside 
#> 225 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 reflect 
#> 226 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 contract inside 
#> 227 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 contract outside 
#> 228 1.10548241135142 0.546223037215453 0.0773191566225074 0.133802817089392 -182.190266239158 contract inside 
#> 229 1.10581486658636 0.548303241564542 0.0780290220884311 0.133891314356184 -182.190262485341 reflect 
#> 230 1.10581486658636 0.548303241564542 0.0780290220884311 0.133891314356184 -182.190262485341 reflect 
#> 231 1.10581486658636 0.548303241564542 0.0780290220884311 0.133891314356184 -182.190262485341 reflect 
#> 232 1.10413999310277 0.545926495456793 0.0776731461298696 0.133982648495521 -182.1902190427 contract inside 
#> 233 1.10413999310277 0.545926495456793 0.0776731461298696 0.133982648495521 -182.1902190427 contract outside 
#> 234 1.10413999310277 0.545926495456793 0.0776731461298696 0.133982648495521 -182.1902190427 contract inside 
#> 235 1.10413999310277 0.545926495456793 0.0776731461298696 0.133982648495521 -182.1902190427 contract inside 
#> 236 1.10559066842464 0.547703666262531 0.0777666680129449 0.133888646841445 -182.190218925431 contract inside 
#> 237 1.10559066842464 0.547703666262531 0.0777666680129449 0.133888646841445 -182.190218925431 reflect 
#> 238 1.10551284045791 0.547257433994207 0.077565763219109 0.133968949876815 -182.190210104138 contract inside 
#> 239 1.10537401554472 0.547241507360136 0.0775420105633108 0.133825411279667 -182.190209868698 contract outside 
#> 240 1.10531101213963 0.54693882618328 0.0775586077956173 0.133883009077873 -182.190208578529 contract inside 
#> 241 1.10479335778349 0.546605604442892 0.077640699371985 0.133937073399036 -182.190207065777 contract inside 
#> 242 1.10479335778349 0.546605604442892 0.077640699371985 0.133937073399036 -182.190207065777 contract outside 
#> 243 1.10479335778349 0.546605604442892 0.077640699371985 0.133937073399036 -182.190207065777 reflect 
#> 244 1.10479335778349 0.546605604442892 0.077640699371985 0.133937073399036 -182.190207065777 reflect 
#> 245 1.10485472401097 0.546522243197586 0.0775544816648026 0.133864174889559 -182.190206524763 contract inside 
#> 246 1.10507087207508 0.546705588141712 0.0775603638182528 0.133898652525002 -182.190206155983 contract inside 
#> 247 1.10495288201847 0.54657362521729 0.0775317956768433 0.133911138175251 -182.19020600784 contract inside 
#> 248 1.10507758891408 0.546854083595769 0.0775720174643955 0.133881714038514 -182.190205879997 contract outside 
#> 249 1.10489118071016 0.54663473884776 0.07759768019741 0.133912996003735 -182.190205564356 contract inside 
#> 250 1.10489118071016 0.54663473884776 0.07759768019741 0.133912996003735 -182.190205564356 reflect 
#> 251 1.10496072305092 0.546756514014763 0.0775786050009735 0.133923309955994 -182.190205372528 reflect 
#> 252 1.10496072305092 0.546756514014763 0.0775786050009735 0.133923309955994 -182.190205372528 contract inside 
#> 253 1.10496072305092 0.546756514014763 0.0775786050009735 0.133923309955994 -182.190205372528 reflect 
#> 254 1.10503939781626 0.54676545399602 0.0775776441060703 0.133932870853573 -182.190205364129 contract inside 
#> 255 1.10503939781626 0.54676545399602 0.0775776441060703 0.133932870853573 -182.190205364129 contract inside 
#> 256 1.10493632836271 0.546674389254703 0.077585444787142 0.133920231244079 -182.190205325091 contract inside 
#> 257 1.10497727443932 0.54669451738025 0.0775684693905302 0.133920980064046 -182.190205279316 contract inside 
#> 258 1.10497727443932 0.54669451738025 0.0775684693905302 0.133920980064046 -182.190205279316 reflect 
#> 259 1.10497653449773 0.546743355140738 0.0775776200865184 0.133921864174865 -182.190205277625 contract inside 
#> 260 1.10500801379174 0.546745062713917 0.0775771366372897 0.133925268759896 -182.190205274301 contract inside 
#> 261 1.10500801379174 0.546745062713917 0.0775771366372897 0.133925268759896 -182.190205274301 contract inside 
#> 262 1.10501581848435 0.54676279974833 0.0775695156283825 0.133920991298002 -182.190205267199 contract outside 
#> 263 1.10499388421046 0.546729460126818 0.0775717405749405 0.133925994725029 -182.190205261187 contract outside 
#> 264 1.10498791848577 0.546719842946333 0.0775712362984654 0.133922254898294 -182.190205260143 contract inside 
#> 265 1.10498897151783 0.546741323173892 0.0775750136755645 0.133922745795011 -182.190205254961 contract inside 
#> 266 1.10498897151783 0.546741323173892 0.0775750136755645 0.133922745795011 -182.190205254961 reflect 
#> 267 1.10498897151783 0.546741323173892 0.0775750136755645 0.133922745795011 -182.190205254961 reflect 
#> 268 1.10498897151783 0.546741323173892 0.0775750136755645 0.133922745795011 -182.190205254961 reflect 
#> 269 1.10498205570312 0.546720850176566 0.0775713062243922 0.133922082622491 -182.190205253417 contract inside 
#> 270 1.10498205570312 0.546720850176566 0.0775713062243922 0.133922082622491 -182.190205253417 reflect 
#> 271 1.10497451303541 0.546714643023501 0.0775719829279634 0.133923867054187 -182.190205253016 contract inside 
#> 272 1.10498682923633 0.546728540136927 0.077571502680116 0.133924131799054 -182.190205252604 contract inside 
#> 273 1.1049841874582 0.546728994532055 0.0775695339272563 0.133921965707865 -182.190205251943 contract inside 
#> 274 1.10498543392711 0.546732290006442 0.0775730475537624 0.133922878795005 -182.190205251899 contract inside 
#> 275 1.10498342610357 0.546731383649312 0.0775717273172266 0.133924339056514 -182.190205251555 reflect 
#> 276 1.10498342610357 0.546731383649312 0.0775717273172266 0.133924339056514 -182.190205251555 reflect 
#> 277 1.10498697372736 0.546731598753244 0.0775714052901779 0.13392356267183 -182.190205251518 contract inside 
#> 278 1.10499021535345 0.546738514057637 0.0775711756646697 0.13392298859048 -182.190205251466 contract inside 
#> 279 1.10498534986593 0.54673122056835 0.0775706864403813 0.133922703992532 -182.190205251223 contract inside 
#> 280 1.10498534986593 0.54673122056835 0.0775706864403813 0.133922703992532 -182.190205251223 contract inside 
#> 281 1.10498527574182 0.546732450322107 0.0775715405972458 0.133923718770504 -182.19020525118 contract inside 
#> 282 1.10498527574182 0.546732450322107 0.0775715405972458 0.133923718770504 -182.19020525118 reflect 
#> 283 1.10498527574182 0.546732450322107 0.0775715405972458 0.133923718770504 -182.19020525118 reflect 
#> 284 1.10498527574182 0.546732450322107 0.0775715405972458 0.133923718770504 -182.19020525118 reflect 
#> 285 1.10498510514931 0.546733196895986 0.0775712379284599 0.133922927935694 -182.190205251127 contract inside 
#> 286 1.10498301406567 0.546729773214402 0.0775713439743343 0.133923118369101 -182.190205251116 contract inside 
#> 287 1.1049847506304 0.546731392956115 0.0775709209011781 0.133922948008462 -182.190205251094 contract inside 
#> 288 1.10498387351107 0.546731272146833 0.0775708798994167 0.133923090646932 -182.190205251094 contract inside 
#> 289 1.10498473079017 0.546731929562069 0.0775713181364822 0.133923370005165 -182.190205251087 contract inside 
#> 290 1.10498358579934 0.546730039508316 0.0775710546274948 0.133923233668273 -182.190205251085 contract outside 
#> 291 1.10498484574163 0.546731851208265 0.0775708930995693 0.13392318168873 -182.190205251077 contract outside 
#> 292 1.10498450479538 0.546731333031056 0.0775709786709419 0.133923083505348 -182.190205251074 contract inside 
#> 293 1.10498450479538 0.546731333031056 0.0775709786709419 0.133923083505348 -182.190205251074 contract inside 
#> 294 1.10498450058033 0.546731527778859 0.0775711461825254 0.133923266601853 -182.190205251072 contract inside 
#> 295 1.10498450058033 0.546731527778859 0.0775711461825254 0.133923266601853 -182.190205251072 contract outside 
#> 296 1.10498450058033 0.546731527778859 0.0775711461825254 0.133923266601853 -182.190205251072 reflect 
#> 297 1.10498450058033 0.546731527778859 0.0775711461825254 0.133923266601853 -182.190205251072 reflect 
#> 298 1.10498457773826 0.546731543820467 0.0775710391837306 0.133923130368397 -182.19020525107 contract inside 
#> 299 1.10498447532125 0.546731570152369 0.0775711067837003 0.133923157594665 -182.19020525107 contract inside 
#> 300 1.10498447532125 0.546731570152369 0.0775711067837003 0.133923157594665 -182.19020525107 contract inside 
#> 301 1.10498447532125 0.546731570152369 0.0775711067837003 0.133923157594665 -182.19020525107 contract outside 
#> 302 1.10498453174208 0.546731582681714 0.0775711026227784 0.133923212639289 -182.190205251069 contract inside 
#> 303 1.10498456457153 0.546731595560175 0.0775710570532634 0.133923154806425 -182.190205251069 contract inside 
#> 304 1.10498449969259 0.546731586038903 0.0775710638682348 0.133923179126548 -182.190205251069 contract inside 
#> 305 1.10498449969259 0.546731586038903 0.0775710638682348 0.133923179126548 -182.190205251069 contract inside 
#> Optimization has terminated successfully.

Now we can see what our results look like.

CR_par_complete <- secsse::extract_par_vals(idparslist, complete_tree_ml_CR$MLpars)
complete_tree_ml_CR
#> $MLpars
#> $MLpars[[1]]
#> $MLpars[[1]][[1]]
#>          0A 1A 0B 1B
#> 0A 1.104984  0  0  0
#> 1A 0.000000  0  0  0
#> 0B 0.000000  0  0  0
#> 1B 0.000000  0  0  0
#> 
#> $MLpars[[1]][[2]]
#>    0A       1A 0B 1B
#> 0A  0 0.000000  0  0
#> 1A  0 1.104984  0  0
#> 0B  0 0.000000  0  0
#> 1B  0 0.000000  0  0
#> 
#> $MLpars[[1]][[3]]
#>    0A 1A       0B 1B
#> 0A  0  0 0.000000  0
#> 1A  0  0 0.000000  0
#> 0B  0  0 1.104984  0
#> 1B  0  0 0.000000  0
#> 
#> $MLpars[[1]][[4]]
#>    0A 1A 0B       1B
#> 0A  0  0  0 0.000000
#> 1A  0  0  0 0.000000
#> 0B  0  0  0 0.000000
#> 1B  0  0  0 1.104984
#> 
#> 
#> $MLpars[[2]]
#>        0A        1A        0B        1B 
#> 0.5467316 0.5467316 0.5467316 0.5467316 
#> 
#> $MLpars[[3]]
#>           0A         1A         0B         1B
#> 0A        NA 0.07757106 0.07757106 0.00000000
#> 1A 0.1339232         NA 0.00000000 0.07757106
#> 0B 0.1339232 0.00000000         NA 0.07757106
#> 1B 0.0000000 0.13392318 0.13392318         NA
#> 
#> 
#> $ML
#> [1] -182.1902
#> 
#> $conv
#> [1] 0
CR_par_complete
#> [1] 1.10498450 0.54673159 0.07757106 0.13392318
spec_rates_complete <- CR_par_complete[1]
ext_rates_complete <- CR_par_complete[2]
Q_01_complete <- CR_par_complete[3]
Q_10_complete <- CR_par_complete[4]
spec_rates_complete
#> [1] 1.104984
ext_rates_complete
#> [1] 0.5467316
Q_01_complete
#> [1] 0.07757106
Q_10_complete
#> [1] 0.1339232

Comparing with reconstructed trees

It would be interesting to see how they compare with the same tree without any extant species. Let’s follow the standard procedure using a similar phylogeny - the same tree we used before - but where all the extinct lineages have been removed. We’ll keep all other model specification the same.


sim_tree_reconstructed <- secsse::secsse_sim(lambdas = sim_lambda_list,
                                             mus = sim_mu_vector,
                                             qs = sim_q_matrix,
                                             crown_age = 5,
                                             num_concealed_states = 2,
                                             seed = 40,
                                             drop_extinct = TRUE)

if (requireNamespace("diversitree")) {
  traits_for_plot_reconstructed <- data.frame(
    trait = as.numeric(sim_tree_reconstructed$obs_traits),
    row.names = sim_tree_reconstructed$phy$tip.label
  )
  diversitree::trait.plot(tree = sim_tree_reconstructed$phy,
                          dat = traits_for_plot_reconstructed,
                          cols = list("trait" = c("blue", "red")),
                          type = "p")
} else {
  plot(sim_tree_reconstructed$phy)
}


reconstructed_tree_ml <- secsse_ml(phy = sim_tree_reconstructed$phy,
                                   traits = sim_tree_reconstructed$obs_traits,
                                   num_concealed_states = 2,
                                   idparslist = idparslist,
                                   idparsopt = idparsopt,
                                   initparsopt = initparsopt,
                                   idparsfix = idparsfix,
                                   parsfix = initparsfix,
                                   sampling_fraction = sampling_fraction,
                                   verbose = FALSE,
                                   is_complete_tree = FALSE)
#> Note: you set some transitions as impossible to happen.
#> 1 0.1 0.1 0.1 0.1 -261.93005389028 initial 
#> 2 0.105 0.1 0.1 0.1 -258.216382647175 reflect 
#> 3 0.105 0.1 0.1 0.1 -258.216382647175 reflect 
#> 4 0.10846394984326 0.0944288126055149 0.0990675332014694 0.10846394984326 -254.998810494622 expand 
#> 5 0.114897886062343 0.0902873160476524 0.0983692199082245 0.101011959521619 -250.4934379125 expand 
#> 6 0.114897886062343 0.0902873160476524 0.0983692199082245 0.101011959521619 -250.4934379125 reflect 
#> 7 0.124172002998582 0.0905356861233287 0.101313961616394 0.102660370839572 -244.749491435056 expand 
#> 8 0.136331330154724 0.0787355101510194 0.0948548055460873 0.112245294610444 -236.655900817446 expand 
#> 9 0.136331330154724 0.0787355101510194 0.0948548055460873 0.112245294610444 -236.655900817446 reflect 
#> 10 0.157463962842126 0.0655976490744996 0.103698336063287 0.104834174336294 -225.674553815757 expand 
#> 11 0.188580384493981 0.0584948324319094 0.0996649403718381 0.11398400869047 -212.855846312177 expand 
#> 12 0.188580384493981 0.0584948324319094 0.0996649403718381 0.11398400869047 -212.855846312177 reflect 
#> 13 0.232811307810817 0.0269459250512547 0.104513689147747 0.130943534064202 -197.615125559306 expand 
#> 14 0.250333207168013 0.0242967293163897 0.106866374342769 0.119320701155624 -193.521807908476 reflect 
#> 15 0.275886868150326 0.0164978330464061 0.0995375741889073 0.134415859843176 -187.693868953716 reflect 
#> 16 0.290567786992314 0.0104768915751455 0.109885889970718 0.13578248786944 -184.819266569516 reflect 
#> 17 0.290567786992314 0.0104768915751455 0.109885889970718 0.13578248786944 -184.819266569516 contract inside 
#> 18 0.302225873216359 0.00107913758124338 0.107970254509982 0.138301143863997 -182.448350614842 reflect 
#> 19 0.302225873216359 0.00107913758124338 0.107970254509982 0.138301143863997 -182.448350614842 contract inside 
#> 20 0.302225873216359 0.00107913758124338 0.107970254509982 0.138301143863997 -182.448350614842 contract inside 
#> 21 0.312285510183205 0.00318560730171753 0.106518346468716 0.136352050445652 -181.019896629056 reflect 
#> 22 0.312285510183205 0.00318560730171753 0.106518346468716 0.136352050445652 -181.019896629056 contract inside 
#> 23 0.312285510183205 0.00318560730171753 0.106518346468716 0.136352050445652 -181.019896629056 contract inside 
#> 24 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -180.372419320842 reflect 
#> 25 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -180.372419320842 contract inside 
#> 26 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -180.372419320842 contract inside 
#> 27 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -180.372419320842 contract inside 
#> 28 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -180.372419320842 contract inside 
#> 29 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -180.372419320842 contract inside 
#> 30 0.315760258880613 0.000340625420629874 0.107848623994427 0.141102405694079 -180.372419320842 contract inside 
#> 31 0.324834680516738 0.00584040227609538 0.105912860281718 0.138329303078879 -179.33715682238 expand 
#> 32 0.334872383886821 0.00197159955397078 0.106621312004214 0.139699044102768 -177.842487507842 expand 
#> 33 0.3549821525049 0.00179668425993489 0.104942650221956 0.140634876279738 -175.394054894606 expand 
#> 34 0.3549821525049 0.00179668425993489 0.104942650221956 0.140634876279738 -175.394054894606 reflect 
#> 35 0.397312367834679 0.00789045650264323 0.102036562199492 0.139448214741894 -171.337954874148 expand 
#> 36 0.396501203659512 0.000878514006177405 0.103954587781011 0.143343896234882 -171.190194491546 reflect 
#> 37 0.463562902846739 0.00530138303480689 0.0995883134492947 0.145849749794543 -166.498946187917 expand 
#> 38 0.463562902846739 0.00530138303480689 0.0995883134492947 0.145849749794543 -166.498946187917 reflect 
#> 39 0.597762149306162 0.0116123084812065 0.0937155563743466 0.146022366563462 -161.99876918793 expand 
#> 40 0.597762149306162 0.0116123084812065 0.0937155563743466 0.146022366563462 -161.99876918793 reflect 
#> 41 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 reflect 
#> 42 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 reflect 
#> 43 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 contract outside 
#> 44 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 contract inside 
#> 45 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 contract inside 
#> 46 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 contract inside 
#> 47 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 contract inside 
#> 48 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 contract inside 
#> 49 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 contract inside 
#> 50 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 contract inside 
#> 51 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 contract inside 
#> 52 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 contract inside 
#> 53 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 contract inside 
#> 54 0.659396562368822 0.0127328085844898 0.0904127751142135 0.147490591928089 -161.507467150133 reflect 
#> 55 0.663423518727895 0.0116529052788619 0.0910290009533167 0.148896745241483 -161.50678667029 contract inside 
#> 56 0.660411161061345 0.0128071163611284 0.0904903931228948 0.147958303901119 -161.506532093247 reflect 
#> 57 0.669183937727868 0.0123983281347956 0.0905099248616501 0.148573019280109 -161.505159682934 contract inside 
#> 58 0.669183937727868 0.0123983281347956 0.0905099248616501 0.148573019280109 -161.505159682934 contract outside 
#> 59 0.669183937727868 0.0123983281347956 0.0905099248616501 0.148573019280109 -161.505159682934 reflect 
#> 60 0.665970017731309 0.0144794412848252 0.0895537368245529 0.146897956392397 -161.499080047311 expand 
#> 61 0.665970017731309 0.0144794412848252 0.0895537368245529 0.146897956392397 -161.499080047311 reflect 
#> 62 0.665970017731309 0.0144794412848252 0.0895537368245529 0.146897956392397 -161.499080047311 contract inside 
#> 63 0.668282055274432 0.0154699305199633 0.0887950645917839 0.145523176083209 -161.496644345179 expand 
#> 64 0.668282055274432 0.0154699305199633 0.0887950645917839 0.145523176083209 -161.496644345179 reflect 
#> 65 0.668282055274432 0.0154699305199633 0.0887950645917839 0.145523176083209 -161.496644345179 reflect 
#> 66 0.667106799215811 0.0203134968489273 0.0863080890862721 0.141591894369165 -161.49258566817 expand 
#> 67 0.667106799215811 0.0203134968489273 0.0863080890862721 0.141591894369165 -161.49258566817 reflect 
#> 68 0.670687354469963 0.0196914814211509 0.0864942042034117 0.141857056648167 -161.492214810743 reflect 
#> 69 0.670687354469963 0.0196914814211509 0.0864942042034117 0.141857056648167 -161.492214810743 reflect 
#> 70 0.670687354469963 0.0196914814211509 0.0864942042034117 0.141857056648167 -161.492214810743 reflect 
#> 71 0.670687354469963 0.0196914814211509 0.0864942042034117 0.141857056648167 -161.492214810743 contract inside 
#> 72 0.667953758520859 0.020578831337367 0.0862399613845515 0.141466138502421 -161.491931283398 reflect 
#> 73 0.670829619420049 0.0199121450085509 0.0864561258547648 0.142012835305563 -161.491688198239 contract inside 
#> 74 0.670829619420049 0.0199121450085509 0.0864561258547648 0.142012835305563 -161.491688198239 contract inside 
#> 75 0.67161620617801 0.0216997654965993 0.0854962223281683 0.140517297788725 -161.49157876483 contract inside 
#> 76 0.668980013132459 0.0216951528457792 0.08566969431932 0.140832086199817 -161.491257845391 reflect 
#> 77 0.668980013132459 0.0216951528457792 0.08566969431932 0.140832086199817 -161.491257845391 reflect 
#> 78 0.668980013132459 0.0216951528457792 0.08566969431932 0.140832086199817 -161.491257845391 reflect 
#> 79 0.668980013132459 0.0216951528457792 0.08566969431932 0.140832086199817 -161.491257845391 contract outside 
#> 80 0.671139803339359 0.0220771430936268 0.0854609478366071 0.140923567432272 -161.490280871432 expand 
#> 81 0.671139803339359 0.0220771430936268 0.0854609478366071 0.140923567432272 -161.490280871432 reflect 
#> 82 0.671139803339359 0.0220771430936268 0.0854609478366071 0.140923567432272 -161.490280871432 contract inside 
#> 83 0.671139803339359 0.0220771430936268 0.0854609478366071 0.140923567432272 -161.490280871432 reflect 
#> 84 0.673711345643402 0.0229308697962357 0.0850914808544844 0.140799640654125 -161.489851839452 expand 
#> 85 0.669785823039254 0.0230447401367683 0.0853007197465296 0.141244709359189 -161.488699252701 expand 
#> 86 0.669785823039254 0.0230447401367683 0.0853007197465296 0.141244709359189 -161.488699252701 reflect 
#> 87 0.669785823039254 0.0230447401367683 0.0853007197465296 0.141244709359189 -161.488699252701 reflect 
#> 88 0.674676572522523 0.0245838769915527 0.0847563549310053 0.141767129899724 -161.486310303238 expand 
#> 89 0.674676572522523 0.0245838769915527 0.0847563549310053 0.141767129899724 -161.486310303238 reflect 
#> 90 0.673775034845731 0.0266175506851315 0.0840062218125982 0.141175869194833 -161.484265326231 expand 
#> 91 0.673775034845731 0.0266175506851315 0.0840062218125982 0.141175869194833 -161.484265326231 reflect 
#> 92 0.678479932314745 0.0292832133766018 0.0833809018115783 0.142855661555204 -161.481005292961 expand 
#> 93 0.678479932314745 0.0292832133766018 0.0833809018115783 0.142855661555204 -161.481005292961 reflect 
#> 94 0.67606807101661 0.0339493993250076 0.0819626073368801 0.14244131696565 -161.47584656177 expand 
#> 95 0.67606807101661 0.0339493993250076 0.0819626073368801 0.14244131696565 -161.47584656177 reflect 
#> 96 0.67606807101661 0.0339493993250076 0.0819626073368801 0.14244131696565 -161.47584656177 reflect 
#> 97 0.685002769877256 0.0419672715728321 0.0795828726420148 0.144408643938863 -161.472766644186 expand 
#> 98 0.687910078302617 0.0517397804661513 0.0764394374702766 0.143554371526049 -161.470928058282 expand 
#> 99 0.684029694161048 0.0481421661636047 0.0779272116835181 0.144850216572704 -161.469875614687 reflect 
#> 100 0.684029694161048 0.0481421661636047 0.0779272116835181 0.144850216572704 -161.469875614687 reflect 
#> 101 0.684029694161048 0.0481421661636047 0.0779272116835181 0.144850216572704 -161.469875614687 contract inside 
#> 102 0.684029694161048 0.0481421661636047 0.0779272116835181 0.144850216572704 -161.469875614687 contract outside 
#> 103 0.684029694161048 0.0481421661636047 0.0779272116835181 0.144850216572704 -161.469875614687 contract inside 
#> 104 0.684029694161048 0.0481421661636047 0.0779272116835181 0.144850216572704 -161.469875614687 contract inside 
#> 105 0.684029694161048 0.0481421661636047 0.0779272116835181 0.144850216572704 -161.469875614687 reflect 
#> 106 0.684029694161048 0.0481421661636047 0.0779272116835181 0.144850216572704 -161.469875614687 contract inside 
#> 107 0.684029694161048 0.0481421661636047 0.0779272116835181 0.144850216572704 -161.469875614687 reflect 
#> 108 0.683837813561719 0.0467179241242624 0.078201026972046 0.144056669037386 -161.469874804306 contract inside 
#> 109 0.684202151519227 0.0486032101989091 0.0778230376781507 0.144875128138363 -161.469741592241 reflect 
#> 110 0.684202151519227 0.0486032101989091 0.0778230376781507 0.144875128138363 -161.469741592241 contract inside 
#> 111 0.684202151519227 0.0486032101989091 0.0778230376781507 0.144875128138363 -161.469741592241 contract inside 
#> 112 0.683339119238023 0.0469035934446589 0.07815744669578 0.143928253068464 -161.469727250192 reflect 
#> 113 0.683281592999339 0.0484214534376091 0.0778618757866948 0.144657336156222 -161.46971576671 reflect 
#> 114 0.683281592999339 0.0484214534376091 0.0778618757866948 0.144657336156222 -161.46971576671 reflect 
#> 115 0.683857464324816 0.0484776587359124 0.0779478697880712 0.144969568922736 -161.469509356635 expand 
#> 116 0.683857464324816 0.0484776587359124 0.0779478697880712 0.144969568922736 -161.469509356635 reflect 
#> 117 0.683857464324816 0.0484776587359124 0.0779478697880712 0.144969568922736 -161.469509356635 contract inside 
#> 118 0.683857464324816 0.0484776587359124 0.0779478697880712 0.144969568922736 -161.469509356635 reflect 
#> 119 0.683304098224695 0.0465173566502098 0.0784163376101304 0.144097874683932 -161.46929750012 expand 
#> 120 0.682694844492124 0.0467516299017768 0.0785686903667074 0.144823636945749 -161.469113498723 expand 
#> 121 0.682694844492124 0.0467516299017768 0.0785686903667074 0.144823636945749 -161.469113498723 reflect 
#> 122 0.684496321225524 0.0475086834871514 0.0784593717927381 0.145148127915802 -161.468805983913 expand 
#> 123 0.684496321225524 0.0475086834871514 0.0784593717927381 0.145148127915802 -161.468805983913 reflect 
#> 124 0.684496321225524 0.0475086834871514 0.0784593717927381 0.145148127915802 -161.468805983913 reflect 
#> 125 0.683991587095499 0.0469405037431861 0.0788279711931932 0.145064421565607 -161.467950866025 expand 
#> 126 0.683991587095499 0.0469405037431861 0.0788279711931932 0.145064421565607 -161.467950866025 reflect 
#> 127 0.683991587095499 0.0469405037431861 0.0788279711931932 0.145064421565607 -161.467950866025 reflect 
#> 128 0.687626740585482 0.0504458204895628 0.078057222380269 0.145836392487511 -161.467913808735 expand 
#> 129 0.686747312028655 0.0474301596946844 0.079221896980219 0.145271129051444 -161.466610956188 expand 
#> 130 0.686747312028655 0.0474301596946844 0.079221896980219 0.145271129051444 -161.466610956188 reflect 
#> 131 0.689473030938786 0.052750665488542 0.0782442989572882 0.146952314390032 -161.466400516365 expand 
#> 132 0.689473030938786 0.052750665488542 0.0782442989572882 0.146952314390032 -161.466400516365 reflect 
#> 133 0.690690152923801 0.0514843341334549 0.079157232601871 0.145902940983361 -161.463545440046 expand 
#> 134 0.690690152923801 0.0514843341334549 0.079157232601871 0.145902940983361 -161.463545440046 reflect 
#> 135 0.690690152923801 0.0514843341334549 0.079157232601871 0.145902940983361 -161.463545440046 reflect 
#> 136 0.690690152923801 0.0514843341334549 0.079157232601871 0.145902940983361 -161.463545440046 reflect 
#> 137 0.694487559696619 0.0509730709400166 0.0813424921243452 0.146933884953498 -161.462905680032 expand 
#> 138 0.690683812050656 0.0497604557531513 0.0812233846117851 0.145625189449406 -161.458770902345 expand 
#> 139 0.690683812050656 0.0497604557531513 0.0812233846117851 0.145625189449406 -161.458770902345 reflect 
#> 140 0.690683812050656 0.0497604557531513 0.0812233846117851 0.145625189449406 -161.458770902345 reflect 
#> 141 0.690683812050656 0.0497604557531513 0.0812233846117851 0.145625189449406 -161.458770902345 reflect 
#> 142 0.690683812050656 0.0497604557531513 0.0812233846117851 0.145625189449406 -161.458770902345 reflect 
#> 143 0.697985997508701 0.0525114743436732 0.084087987067195 0.145386628931684 -161.456381330457 expand 
#> 144 0.692229992041617 0.0547843298602958 0.0797764509783189 0.14421939712557 -161.45491644665 expand 
#> 145 0.689110735154414 0.0516718102788699 0.0823166462094521 0.143714827731044 -161.449682449329 expand 
#> 146 0.700165131106412 0.0679503826808044 0.0792635003322385 0.147998605025567 -161.448613907665 expand 
#> 147 0.703279039129958 0.0708235885636648 0.0816261904164096 0.144732505363657 -161.433811439507 expand 
#> 148 0.703279039129958 0.0708235885636648 0.0816261904164096 0.144732505363657 -161.433811439507 reflect 
#> 149 0.705754115105173 0.0862980194414938 0.0809079108274218 0.147601542693936 -161.4218272534 expand 
#> 150 0.705754115105173 0.0862980194414938 0.0809079108274218 0.147601542693936 -161.4218272534 reflect 
#> 151 0.711812996082886 0.107731691001953 0.0794031502884629 0.143659514613303 -161.410072399827 expand 
#> 152 0.73718629767034 0.13282132774301 0.0846096959389692 0.148807171383058 -161.399298051218 expand 
#> 153 0.715991516433322 0.101087669083292 0.0860024035813745 0.143475661556152 -161.396954698638 reflect 
#> 154 0.732174376430318 0.145129906878294 0.0838243996929365 0.147032043251543 -161.392908609684 reflect 
#> 155 0.743101567072381 0.158857118877991 0.0860126693606941 0.143882731426102 -161.387956899728 reflect 
#> 156 0.743101567072381 0.158857118877991 0.0860126693606941 0.143882731426102 -161.387956899728 reflect 
#> 157 0.734623936976011 0.149630219005442 0.0887459024506324 0.142368956964346 -161.385490239495 reflect 
#> 158 0.734623936976011 0.149630219005442 0.0887459024506324 0.142368956964346 -161.385490239495 contract inside 
#> 159 0.734623936976011 0.149630219005442 0.0887459024506324 0.142368956964346 -161.385490239495 reflect 
#> 160 0.734623936976011 0.149630219005442 0.0887459024506324 0.142368956964346 -161.385490239495 contract inside 
#> 161 0.74248181797864 0.1501567706403 0.090062777949647 0.143155531424965 -161.384978047288 contract inside 
#> 162 0.732285907672949 0.131703011978217 0.0915735609041636 0.143859000094254 -161.383950975406 reflect 
#> 163 0.739917961354324 0.146831981769975 0.089461931896396 0.144508400194276 -161.383391274475 contract inside 
#> 164 0.732759151066328 0.135606872247579 0.088317764942564 0.143929857138088 -161.38317720668 contract inside 
#> 165 0.735736197690808 0.145310702174855 0.0892990422190331 0.143115488203076 -161.38297761782 contract inside 
#> 166 0.727918258413732 0.129681873385312 0.0892610841799782 0.144551264761694 -161.382599706009 reflect 
#> 167 0.727918258413732 0.129681873385312 0.0892610841799782 0.144551264761694 -161.382599706009 contract inside 
#> 168 0.736147523700619 0.141640690745732 0.0893814325136909 0.144196386638374 -161.382595868716 contract inside 
#> 169 0.733718889495823 0.14040526220433 0.0908194019635511 0.143972460880038 -161.382011298175 reflect 
#> 170 0.733718889495823 0.14040526220433 0.0908194019635511 0.143972460880038 -161.382011298175 contract inside 
#> 171 0.733718889495823 0.14040526220433 0.0908194019635511 0.143972460880038 -161.382011298175 reflect 
#> 172 0.733718889495823 0.14040526220433 0.0908194019635511 0.143972460880038 -161.382011298175 contract inside 
#> 173 0.733718889495823 0.14040526220433 0.0908194019635511 0.143972460880038 -161.382011298175 reflect 
#> 174 0.729464684859982 0.135829962508872 0.0902493924299325 0.144540417616807 -161.381882271972 reflect 
#> 175 0.729464684859982 0.135829962508872 0.0902493924299325 0.144540417616807 -161.381882271972 reflect 
#> 176 0.729464684859982 0.135829962508872 0.0902493924299325 0.144540417616807 -161.381882271972 reflect 
#> 177 0.732987680518715 0.139547002440537 0.09133408458326 0.144539725182982 -161.381784383964 reflect 
#> 178 0.732987680518715 0.139547002440537 0.09133408458326 0.144539725182982 -161.381784383964 reflect 
#> 179 0.731273132658095 0.138285312238593 0.0898819296371555 0.145074839057771 -161.381755791056 reflect 
#> 180 0.732814081335451 0.139094207863766 0.0900457967605005 0.144638535271112 -161.381736173161 contract inside 
#> 181 0.732814081335451 0.139094207863766 0.0900457967605005 0.144638535271112 -161.381736173161 reflect 
#> 182 0.732814081335451 0.139094207863766 0.0900457967605005 0.144638535271112 -161.381736173161 reflect 
#> 183 0.731105930899161 0.136237451723976 0.0905594825255336 0.144934327563819 -161.381657569373 contract inside 
#> 184 0.732620237764866 0.138892026817532 0.0908146344480089 0.144738431527486 -161.381650268911 contract inside 
#> 185 0.731931271930235 0.138337358268712 0.0902053099327837 0.144963962134699 -161.381645859617 contract inside 
#> 186 0.731931271930235 0.138337358268712 0.0902053099327837 0.144963962134699 -161.381645859617 contract inside 
#> 187 0.731499141331344 0.137008613692747 0.0910197943689705 0.14515993050894 -161.381612507997 reflect 
#> 188 0.731499141331344 0.137008613692747 0.0910197943689705 0.14515993050894 -161.381612507997 reflect 
#> 189 0.732224319945746 0.139133086614274 0.0908344560067091 0.14496593784696 -161.381600326734 reflect 
#> 190 0.732224319945746 0.139133086614274 0.0908344560067091 0.14496593784696 -161.381600326734 reflect 
#> 191 0.732224319945746 0.139133086614274 0.0908344560067091 0.14496593784696 -161.381600326734 contract inside 
#> 192 0.7322939886238 0.138782544967112 0.090724020779115 0.145274771348808 -161.38158261422 reflect 
#> 193 0.7322939886238 0.138782544967112 0.090724020779115 0.145274771348808 -161.38158261422 contract outside 
#> 194 0.7322939886238 0.138782544967112 0.090724020779115 0.145274771348808 -161.38158261422 contract inside 
#> 195 0.7322939886238 0.138782544967112 0.090724020779115 0.145274771348808 -161.38158261422 contract outside 
#> 196 0.7322939886238 0.138782544967112 0.090724020779115 0.145274771348808 -161.38158261422 contract inside 
#> 197 0.731894784579091 0.13834263379254 0.0908551998217158 0.145263849759435 -161.381579239644 reflect 
#> 198 0.731894784579091 0.13834263379254 0.0908551998217158 0.145263849759435 -161.381579239644 reflect 
#> 199 0.731894784579091 0.13834263379254 0.0908551998217158 0.145263849759435 -161.381579239644 reflect 
#> 200 0.731894784579091 0.13834263379254 0.0908551998217158 0.145263849759435 -161.381579239644 reflect 
#> 201 0.731894784579091 0.13834263379254 0.0908551998217158 0.145263849759435 -161.381579239644 contract inside 
#> 202 0.731894784579091 0.13834263379254 0.0908551998217158 0.145263849759435 -161.381579239644 contract inside 
#> 203 0.731894784579091 0.13834263379254 0.0908551998217158 0.145263849759435 -161.381579239644 contract inside 
#> 204 0.731894784579091 0.13834263379254 0.0908551998217158 0.145263849759435 -161.381579239644 contract inside 
#> 205 0.731894784579091 0.13834263379254 0.0908551998217158 0.145263849759435 -161.381579239644 reflect 
#> 206 0.73212565104751 0.138612684657158 0.090852058615267 0.14531708069411 -161.381579073491 reflect 
#> 207 0.73229107909569 0.138934114766673 0.0908392291904988 0.145241003097985 -161.381578654049 contract inside 
#> 208 0.73229107909569 0.138934114766673 0.0908392291904988 0.145241003097985 -161.381578654049 reflect 
#> 209 0.732141338548979 0.138627403468732 0.0908838390492094 0.145234429220726 -161.381578403872 contract inside 
#> 210 0.732056908034154 0.138537421805048 0.0909037017092758 0.145262967184396 -161.381578390252 contract inside 
#> 211 0.732024252595435 0.138510235185199 0.0908624531357244 0.145263859443586 -161.381578316272 contract inside 
#> 212 0.732127019733565 0.138632476438493 0.0908621818345963 0.145283821675131 -161.381578262441 contract inside 
#> 213 0.732189222745081 0.13875547020094 0.0908586360801145 0.14525113601463 -161.381578110648 contract inside 
#> 214 0.732189222745081 0.13875547020094 0.0908586360801145 0.14525113601463 -161.381578110648 contract inside 
#> 215 0.732189222745081 0.13875547020094 0.0908586360801145 0.14525113601463 -161.381578110648 contract outside 
#> 216 0.732189222745081 0.13875547020094 0.0908586360801145 0.14525113601463 -161.381578110648 contract inside 
#> 217 0.732138488969956 0.138673293045704 0.0908605154048685 0.145242704201734 -161.381578104661 contract outside 
#> 218 0.73212977273927 0.138645810327614 0.0908672714035133 0.145252269372515 -161.381578084088 contract inside 
#> 219 0.732131154374418 0.138661844808451 0.0908700628267886 0.145247428940358 -161.381578083436 contract outside 
#> 220 0.732115950138695 0.138637167627197 0.0908629650306962 0.145255575471004 -161.381578068481 contract inside 
#> 221 0.732115950138695 0.138637167627197 0.0908629650306962 0.145255575471004 -161.381578068481 contract inside 
#> 222 0.732129465199817 0.138651615971764 0.0908705941809536 0.145260090429113 -161.381578063375 reflect 
#> 223 0.732142155144132 0.138700098630744 0.0908646135757272 0.1452555187629 -161.381578043351 expand 
#> 224 0.732142155144132 0.138700098630744 0.0908646135757272 0.1452555187629 -161.381578043351 reflect 
#> 225 0.732142155144132 0.138700098630744 0.0908646135757272 0.1452555187629 -161.381578043351 contract inside 
#> 226 0.732142155144132 0.138700098630744 0.0908646135757272 0.1452555187629 -161.381578043351 reflect 
#> 227 0.732142155144132 0.138700098630744 0.0908646135757272 0.1452555187629 -161.381578043351 contract inside 
#> 228 0.732142155144132 0.138700098630744 0.0908646135757272 0.1452555187629 -161.381578043351 contract outside 
#> 229 0.732142155144132 0.138700098630744 0.0908646135757272 0.1452555187629 -161.381578043351 contract inside 
#> 230 0.732147066216019 0.13870619819624 0.0908695626893005 0.145259632892027 -161.381578040843 reflect 
#> 231 0.732147066216019 0.13870619819624 0.0908695626893005 0.145259632892027 -161.381578040843 contract inside 
#> 232 0.732147066216019 0.13870619819624 0.0908695626893005 0.145259632892027 -161.381578040843 reflect 
#> 233 0.732147066216019 0.13870619819624 0.0908695626893005 0.145259632892027 -161.381578040843 reflect 
#> 234 0.732142047477659 0.138698214272073 0.090866021937515 0.14525748620451 -161.38157804039 contract inside 
#> 235 0.732142047477659 0.138698214272073 0.090866021937515 0.14525748620451 -161.38157804039 reflect 
#> 236 0.732142265147751 0.138693244210182 0.0908675731673758 0.145258076832295 -161.381578039551 contract inside 
#> 237 0.732142265147751 0.138693244210182 0.0908675731673758 0.145258076832295 -161.381578039551 contract inside 
#> 238 0.732138490825861 0.138693984609258 0.0908676896260855 0.145258407269734 -161.381578039538 contract inside 
#> 239 0.732138490825861 0.138693984609258 0.0908676896260855 0.145258407269734 -161.381578039538 contract inside 
#> 240 0.732140730924374 0.138693070952303 0.0908701222363819 0.145258685280228 -161.381578039514 reflect 
#> 241 0.732141920789993 0.138695737779923 0.0908684271428299 0.145259785949195 -161.381578039278 reflect 
#> 242 0.732137718396679 0.13868739915272 0.0908683993048139 0.145258770916654 -161.381578039085 reflect 
#> 243 0.732137718396679 0.13868739915272 0.0908683993048139 0.145258770916654 -161.381578039085 contract inside 
#> 244 0.732139415449253 0.138693130305878 0.0908682279444694 0.145258670727063 -161.38157803905 contract inside 
#> 245 0.732139415449253 0.138693130305878 0.0908682279444694 0.145258670727063 -161.38157803905 contract inside 
#> 246 0.732140772281255 0.138693632198334 0.0908684574568309 0.14525923599305 -161.381578039023 contract inside 
#> 247 0.732140772281255 0.138693632198334 0.0908684574568309 0.14525923599305 -161.381578039023 reflect 
#> 248 0.732140772281255 0.138693632198334 0.0908684574568309 0.14525923599305 -161.381578039023 reflect 
#> 249 0.73213835824032 0.138689568331948 0.0908683952599928 0.145258924195615 -161.381578038964 contract inside 
#> 250 0.73213835824032 0.138689568331948 0.0908683952599928 0.145258924195615 -161.381578038964 contract inside 
#> 251 0.73213835824032 0.138689568331948 0.0908683952599928 0.145258924195615 -161.381578038964 contract inside 
#> 252 0.73213835824032 0.138689568331948 0.0908683952599928 0.145258924195615 -161.381578038964 contract inside 
#> 253 0.73213835824032 0.138689568331948 0.0908683952599928 0.145258924195615 -161.381578038964 contract inside 
#> 254 0.73213835824032 0.138689568331948 0.0908683952599928 0.145258924195615 -161.381578038964 reflect 
#> 255 0.732139118525708 0.138691633659878 0.0908684714645673 0.145258967099103 -161.381578038956 contract inside 
#> 256 0.732139118525708 0.138691633659878 0.0908684714645673 0.145258967099103 -161.381578038956 reflect 
#> 257 0.732139379315965 0.138691520266516 0.0908685842898004 0.145258960793962 -161.381578038953 contract inside 
#> 258 0.732139425354361 0.138691659726148 0.0908684410696342 0.145259024200888 -161.381578038952 contract inside 
#> 259 0.732138881909154 0.138690585598648 0.0908684270349462 0.145258936606202 -161.381578038949 contract inside 
#> 260 0.732138881909154 0.138690585598648 0.0908684270349462 0.145258936606202 -161.381578038949 contract inside 
#> 261 0.732138881909154 0.138690585598648 0.0908684270349462 0.145258936606202 -161.381578038949 contract inside 
#> 262 0.732138881909154 0.138690585598648 0.0908684270349462 0.145258936606202 -161.381578038949 reflect 
#> 263 0.732138881909154 0.138690585598648 0.0908684270349462 0.145258936606202 -161.381578038949 reflect 
#> 264 0.732139234446918 0.138691208904998 0.0908683973886132 0.145258917824323 -161.381578038948 contract inside 
#> 265 0.732138997939025 0.138690864404218 0.0908683757044896 0.145258899761607 -161.381578038948 contract inside 
#> 266 0.732139089503891 0.138691053923185 0.09086835300463 0.145258941944443 -161.381578038948 contract inside 
#> 267 0.732139126154796 0.138691199981566 0.0908684283880087 0.145258943146545 -161.381578038948 contract inside 
#> 268 0.732139126154796 0.138691199981566 0.0908684283880087 0.145258943146545 -161.381578038948 contract inside 
#> 269 0.732139126154796 0.138691199981566 0.0908684283880087 0.145258943146545 -161.381578038948 reflect 
#> 270 0.732139126154796 0.138691199981566 0.0908684283880087 0.145258943146545 -161.381578038948 reflect 
#> 271 0.732139126154796 0.138691199981566 0.0908684283880087 0.145258943146545 -161.381578038948 contract inside 
#> 272 0.732139126154796 0.138691199981566 0.0908684283880087 0.145258943146545 -161.381578038948 reflect 
#> 273 0.73213901286607 0.138690930710058 0.0908684163623314 0.145258943397563 -161.381578038947 contract inside 
#> 274 0.73213901286607 0.138690930710058 0.0908684163623314 0.145258943397563 -161.381578038947 contract inside 
#> 275 0.732139020898734 0.138690975219647 0.0908684231683261 0.145258930325078 -161.381578038947 contract outside 
#> 276 0.732139020898734 0.138690975219647 0.0908684231683261 0.145258930325078 -161.381578038947 contract inside 
#> 277 0.732139068835574 0.13869108164598 0.0908684253137868 0.145258943053262 -161.381578038947 contract inside 
#> 278 0.732139068835574 0.13869108164598 0.0908684253137868 0.145258943053262 -161.381578038947 contract inside 
#> 279 0.732139068835574 0.13869108164598 0.0908684253137868 0.145258943053262 -161.381578038947 reflect 
#> 280 0.732139068835574 0.13869108164598 0.0908684253137868 0.145258943053262 -161.381578038947 contract inside 
#> Optimization has terminated successfully.
reconstructed_tree_ml_CR <- reconstructed_tree_ml$ML
CR_par_reconstructed <- secsse::extract_par_vals(
  idparslist,
  reconstructed_tree_ml$MLpars
)
reconstructed_tree_ml
#> $MLpars
#> $MLpars[[1]]
#> $MLpars[[1]][[1]]
#>           0A 1A 0B 1B
#> 0A 0.7321391  0  0  0
#> 1A 0.0000000  0  0  0
#> 0B 0.0000000  0  0  0
#> 1B 0.0000000  0  0  0
#> 
#> $MLpars[[1]][[2]]
#>    0A        1A 0B 1B
#> 0A  0 0.0000000  0  0
#> 1A  0 0.7321391  0  0
#> 0B  0 0.0000000  0  0
#> 1B  0 0.0000000  0  0
#> 
#> $MLpars[[1]][[3]]
#>    0A 1A        0B 1B
#> 0A  0  0 0.0000000  0
#> 1A  0  0 0.0000000  0
#> 0B  0  0 0.7321391  0
#> 1B  0  0 0.0000000  0
#> 
#> $MLpars[[1]][[4]]
#>    0A 1A 0B        1B
#> 0A  0  0  0 0.0000000
#> 1A  0  0  0 0.0000000
#> 0B  0  0  0 0.0000000
#> 1B  0  0  0 0.7321391
#> 
#> 
#> $MLpars[[2]]
#>        0A        1A        0B        1B 
#> 0.1386911 0.1386911 0.1386911 0.1386911 
#> 
#> $MLpars[[3]]
#>           0A         1A         0B         1B
#> 0A        NA 0.09086843 0.09086843 0.00000000
#> 1A 0.1452589         NA 0.00000000 0.09086843
#> 0B 0.1452589 0.00000000         NA 0.09086843
#> 1B 0.0000000 0.14525894 0.14525894         NA
#> 
#> 
#> $ML
#> [1] -161.3816
#> 
#> $conv
#> [1] 0
CR_par_reconstructed
#> [1] 0.73213907 0.13869108 0.09086843 0.14525894
spec_rates_reconstructed <- CR_par_reconstructed[1]
ext_rates_reconstructed <- CR_par_reconstructed[2]
Q_01_reconstructed <- CR_par_reconstructed[3]
Q_10_reconstructed <- CR_par_reconstructed[4]

knitr::kable(
  data.frame(
    Reconstructed_tree = c(
      spec_rates_reconstructed,
      ext_rates_reconstructed,
      Q_01_reconstructed,
      Q_10_reconstructed
    ),  
    Complete_tree = c(
      spec_rates_complete,
      ext_rates_complete,
      Q_01_complete,
      Q_10_complete
    ),
    Generating_parameters = c(
      speciation_rate,
      extinction_rate,
      q_01,
      q_10
    ),
    row.names = c(
      "Speciation rate",
      "Extinction rate",
      "Transition rate 01",
      "Transition rate 10"
    )
  )
)
Reconstructed_tree Complete_tree Generating_parameters
Speciation rate 0.7321391 1.1049845 0.8
Extinction rate 0.1386911 0.5467316 0.2
Transition rate 01 0.0908684 0.0775711 0.1
Transition rate 10 0.1452589 0.1339232 0.1

We see that including extinct species results in a better esimation particularly of the extinction rate. This effect is especially noticeable if there are many extinct species present in the tree. Additionally, we see that the estimation of the transition rate from state 1 to 0 also improved.

As a final note, do note that this is just a simple simulation example and care should be exercised with model selection and specification when fitting secsse to empirical datasets to make predictions about evolutionary patterns.

References

Nee S, May RM, Harvey PH. The reconstructed evolutionary process. Philos Trans R Soc Lond B Biol Sci. 1994 May 28;344(1309):305-11. https://doi.org/10.1098/rstb.1994.0068.