### define two models...
birth.logis <- function(a, b) 1/(1 + exp(-a*t - b)) # logistic
birth.step <- function(l1, l2, Tcl) { # 2 rates with one break-point
ans <- rep(l1, length(t))
ans[t > Tcl] <- l2
ans
}
### ... and their primitives:
BIRTH.logis <- function(t) log(exp(-a*t) + exp(b))/a + t
BIRTH.step <- function(t)
{
if (t <= Tcl) return(t*l1)
Tcl*l1 + (t - Tcl)*l2
}
data(bird.families)
### fit both models:
yule.time(bird.families, birth.logis)
yule.time(bird.families, birth.logis, BIRTH.logis) # same but faster
yule.time(bird.families, birth.step) # fails
yule.time(bird.families, birth.step, BIRTH.step)
Run the code above in your browser using DataCamp Workspace