Perform optimization to find a maximum likelihood estimate.

Description

Perform optimization to find a maximum likelihood estimate.

Usage

expoTree.optim(forest,fix=rep(5,F),fix.val=rep(5,0),
      pars=vector(length=sum(!fix))*1,lo=rep(sum(!fix),0),hi=rep(sum(!fix),0),
      survival=TRUE,vflag=0,method="Nelder-Mead",control=list(trace=0)) 

Arguments

Details

Special case: It is possible to fix the ratio between mu and psi. If mu is fixed at a negtive value, it is seen as a fixed ratio. Example: r <- 0.5 fix <- c(F,F,T,F,T) fix.val <- c(0,0,-r,0,0) At every evaluation, mu is calculated accordingly: psi/(psi+mu) = r ==> mu = psi*(1/r-1)

Value

Output from the optimization. List with the following entries: par : parameter estimates value : negative log-likelihood For the other returned values, see help(optim) respectively

Author(s)

Gabriel E. Leventhal

References

Leventhal, Guenthard, Bonhoeffer & Stadler, 2013

See Also

expoTree

Examples

  # simulate trees
  N <- 15
  beta <- 1
  mu <- 0.1
  psi <- 0.1
  rho <- 0
  nsamp <- 20
  epi <- sim.epi(N,beta,mu,psi,nsamp)
  tree <- cbind(epi$times,epi$ttypes)

  extant <- sum(2*tree[,2]-1)
  lineages <- sum(2*tree[,2]-1)+cumsum(1-2*tree[,2])
  max.lineages <- max(lineages)

  # calculate likelihood for the forest
  lik <- runExpoTree(pars=c(N,beta,mu,psi,rho),times=tree[,1],ttypes=tree[,2],
                     survival=TRUE,return.full=FALSE)
  cat("Likelihood = ",lik,"\n")

  if (! any(is.nan(lik))) {
    expoTree.optim(list(tree),pars=c(N,beta,psi),
                   lo=c(max(max.lineages),0,0),hi=c(50,3,1),
                   fix=c(FALSE,FALSE,TRUE,FALSE,TRUE),
                   fix.val=c(0,0,-psi/(psi+mu),0,0))
  }

Results