Skip to contents

Estimate the upper age distribution and output trees for MCMCtree input

Usage

estimateUpper(
  maxAge,
  rightTail = 0.025,
  phy,
  monoGroups,
  writeMCMCtree = FALSE,
  MCMCtreeName = "estimateUpper.tre"
)

Arguments

maxAge

vector of maximum age bounds for nodes matching order in monoGroups

rightTail

probability of right tail (maximum bound default = 0.025)

phy

fully resolved phylogeny in ape format

monoGroups

list with each element containing species that define a node of interest

writeMCMCtree

logical whether to write tree in format that is compatible with MCMCTree to file

MCMCtreeName

MCMCtree.output file name

Value

list containing node estimates for each distribution

  • "parameters" estimated parameters for each node

  • "apePhy" phylogeny in APE format with node labels showing node distributions

  • "MCMCtree" phylogeny in MCMCtree format

  • "nodeLabels" node labels in MCMCtreeR format

If plot=TRUE plot of distributions in file 'pdfOutput' written to current working directory

If writeMCMCtree=TRUE tree in MCMCtree format in file "MCMCtreeName" written to current working directory

Author

Mark Puttick

Examples

data(apeData)
attach(apeData)
#> The following objects are masked from apeData (pos = 3):
#> 
#>     apeTree, maximumTimes, minimumTimes, monophyleticGroups
#> The following objects are masked from apeData (pos = 4):
#> 
#>     apeTree, maximumTimes, minimumTimes, monophyleticGroups
#> The following objects are masked from apeData (pos = 5):
#> 
#>     apeTree, maximumTimes, minimumTimes, monophyleticGroups
#> The following objects are masked from apeData (pos = 6):
#> 
#>     apeTree, maximumTimes, minimumTimes, monophyleticGroups
#> The following objects are masked from apeData (pos = 7):
#> 
#>     apeTree, maximumTimes, minimumTimes, monophyleticGroups
#> The following objects are masked from apeData (pos = 8):
#> 
#>     apeTree, maximumTimes, minimumTimes, monophyleticGroups
## extract taxon descending from calibrated nodes 8, 10, 11, 13
## these nodes can be visualised using plot.phylo
## and nodelabels from ape
monophyleticGroups <- tipDes(apeData$apeTree, c(8,10,11,13))
minimumTimes <- c("nodeOne"=15, "nodeTwo"=6,
"nodeThree"=8, "nodeFour"=13) / 10
maximumTimes <- c("nodeOne" = 30, "nodeTwo" = 12,
"nodeThree"=12, "nodeFour" = 20) / 10
estimateUpper(maxAge=maximumTimes, monoGroups=monophyleticGroups,
rightTail=0.025, phy=apeTree)
#> [1] "warning - maxProb parameter value recycled"
#> Warning: Length of node.label does not match number of nodes.
#> $apePhy
#> 
#> Phylogenetic tree with 7 tips and 6 internal nodes.
#> 
#> Tip labels:
#>   human, chimpanzee, bonobo, gorilla, orangutan, sumatran, ...
#> Node labels:
#>   'U[3~0.025]', NA, 'U[1.2~0.025]', 'U[1.2~0.025]', NA, 'U[2~0.025]'
#> 
#> Rooted; no branch lengths.
#> 
#> $MCMCtree
#>      
#> 1 7 1
#>                                                                                                                             
#> 1 ((((human,(chimpanzee,bonobo))'U(1.2,0.025)',gorilla)'U(1.2,0.025)',(orangutan,sumatran)'U(2,0.025)'),gibbon)'U(3,0.025)';
#>                
#> 1 //end of file
#> 
#> $parameters
#>         tU    pR
#> node_1 3.0 0.025
#> node_2 1.2 0.025
#> node_3 1.2 0.025
#> node_4 2.0 0.025
#> 
#> $nodeLabels
#> [1] "'U[3~0.025]'"   "'U[1.2~0.025]'" "'U[1.2~0.025]'" "'U[2~0.025]'"  
#>