Monday, August 17, 2015

More updates associated with mapped discrete state ("simmap") object class

In the past few days I have been busy doing things like adding S3 methods for the object class "simmap" that I have now created for the myriad of functions in phytools that create or use modified "phylo" objects of this type. I have also posted a new version of phytools with these methods & other updates & bug-fixes associated with the new object class for a variety of functions in the phytools package. I'm hoping to get this version on CRAN & relatively bug-free in time for my workshop next week at Universidad Nacional Autónoma de Mexico (UNAM) in Mexico City next week.

Here is a little bit of a preview:

install.packages("http://www.phytools.org/nonstatic/phytools_0.4-99.tar.gz",
type="source",repos=NULL)
## Installing package into 'C:/Users/Liam/Documents/R/win-library/3.2'
## (as 'lib' is unspecified)
library(phytools)
packageVersion("phytools")
## [1] '0.4.99'
## load "simmap" style tree from package data
data(anoletree)
## S3 print method
anoletree
## 
## Phylogenetic tree with 82 tips and 81 internal nodes.
## 
## Tip labels:
##  Anolis_ahli, Anolis_allogus, Anolis_rubribarbus, Anolis_imias, Anolis_sagrei, Anolis_bremeri, ...
## 
## The tree includes a mapped, 6-state discrete character with states:
##  CG, GB, TC, TG, Tr, Tw
## 
## Rooted; includes branch lengths.
## summary method for object of class "simmap"
summary(anoletree)
## 1 tree with a mapped discrete character with states:
##  CG, GB, TC, TG, Tr, Tw 
## 
## tree has 24 changes between states
## 
## changes are of the following types:
##    CG GB TC TG Tr Tw
## CG  0  0  0  0  0  1
## GB  1  0  0  0  0  0
## TC  1  0  0  0  0  2
## TG  2  6  4  0  1  3
## Tr  0  0  0  0  0  0
## Tw  0  0  1  1  1  0
## 
## mean total time spent in each state is:
##               CG         GB         TC         TG          Tr         Tw
## raw  14.49251700 40.1094769 31.1538539 80.9012968 11.83957979 27.1702360
## prop  0.07046595  0.1950215  0.1514772  0.3933607  0.05756676  0.1321079
##        total
## raw  205.667
## prop   1.000
## plot method
plot(anoletree,fsize=0.6,ftype="i",xlim=c(-0.18*max(nodeHeights(anoletree)),
    1.22*max(nodeHeights(anoletree))))
## no colors provided. using the following legend:
##        CG        GB        TC        TG        Tr        Tw 
##   "black"     "red"  "green3"    "blue"    "cyan" "magenta"
states<-sort(unique(getStates(anoletree,"tips")))
add.simmap.legend(x=-0.18*max(nodeHeights(anoletree)),y=Ntip(anoletree),
    colors=setNames(palette()[1:length(states)],states),prompt=FALSE)

plot of chunk unnamed-chunk-1

## do stochastic character mapping
x<-getStates(anoletree,"tips")
x
##            Anolis_ahli         Anolis_allogus     Anolis_rubribarbus 
##                   "TG"                   "TG"                   "TG" 
##           Anolis_imias          Anolis_sagrei         Anolis_bremeri 
##                   "TG"                   "TG"                   "TG" 
## Anolis_quadriocellifer      Anolis_ophiolepis         Anolis_mestrei 
##                   "TG"                   "GB"                   "TG" 
##
##                 ......
##
##      Anolis_singularis    Anolis_chlorocyanus     Anolis_coelestinus 
##                   "TC"                   "TC"                   "TC" 
##        Anolis_occultus 
##                   "Tw"
maps<-make.simmap(anoletree,x,nsim=100,model="ER")
## make.simmap is sampling character histories conditioned on the transition matrix
## Q =
##             CG          GB          TC          TG          Tr          Tw
## CG -0.11570723  0.02314145  0.02314145  0.02314145  0.02314145  0.02314145
## GB  0.02314145 -0.11570723  0.02314145  0.02314145  0.02314145  0.02314145
## TC  0.02314145  0.02314145 -0.11570723  0.02314145  0.02314145  0.02314145
## TG  0.02314145  0.02314145  0.02314145 -0.11570723  0.02314145  0.02314145
## Tr  0.02314145  0.02314145  0.02314145  0.02314145 -0.11570723  0.02314145
## Tw  0.02314145  0.02314145  0.02314145  0.02314145  0.02314145 -0.11570723
## (estimated using likelihood);
## and (mean) root node prior probabilities
## pi =
##        CG        GB        TC        TG        Tr        Tw 
## 0.1666667 0.1666667 0.1666667 0.1666667 0.1666667 0.1666667
## Done.
maps
## 100 phylogenetic trees with mapped discrete characters
print(maps,details=TRUE)
## 100 phylogenetic trees with mapped discrete characters
## tree 1 : 82 tips, 6 mapped states
## tree 2 : 82 tips, 6 mapped states
## tree 3 : 82 tips, 6 mapped states
## tree 4 : 82 tips, 6 mapped states
## tree 5 : 82 tips, 6 mapped states
## tree 6 : 82 tips, 6 mapped states
## ....
## tree 97 : 82 tips, 6 mapped states
## tree 98 : 82 tips, 6 mapped states
## tree 99 : 82 tips, 6 mapped states
## tree 100 : 82 tips, 6 mapped states
obj<-summary(maps)
plot(obj,fsize=0.6,ftype="i",xlim=c(-0.18*max(nodeHeights(anoletree)),
    1.22*max(nodeHeights(anoletree))))
add.simmap.legend(x=-0.18*max(nodeHeights(anoletree)),y=Ntip(anoletree),
    colors=setNames(palette()[1:length(states)],states),prompt=FALSE)

plot of chunk unnamed-chunk-2

Or, alternatively:

plot(anoletree,fsize=0.6,ftype="i",xlim=c(-0.18*max(nodeHeights(anoletree)),
    1.22*max(nodeHeights(anoletree))))
## no colors provided. using the following legend:
##        CG        GB        TC        TG        Tr        Tw 
##   "black"     "red"  "green3"    "blue"    "cyan" "magenta"
states<-sort(unique(getStates(anoletree,"tips")))
add.simmap.legend(x=-0.18*max(nodeHeights(anoletree)),y=Ntip(anoletree),
    colors=setNames(palette()[1:length(states)],states),prompt=FALSE)
nodelabels(pie=obj$ace,piecol=setNames(palette()[1:length(states)],states),
    cex=0.6)

plot of chunk unnamed-chunk-3

The same methods also apply to "simmap" class objects created by other methods. For instance:

tree<-pbtree(n=50,scale=100)
limits=c(0,25,50,75)
tree<-make.era.map(tree,limits)
tree
## 
## Phylogenetic tree with 50 tips and 49 internal nodes.
## 
## Tip labels:
##  t1, t33, t34, t20, t21, t16, ...
## 
## The tree includes a mapped, 4-state discrete character with states:
##  1, 2, 3, 4
## 
## Rooted; includes branch lengths.
plot(tree,fsize=0.7,mar=c(5.1,0.1,0.1,0.1))
## no colors provided. using the following legend:
##        1        2        3        4 
##  "black"    "red" "green3"   "blue"
axis(1)
title(xlab="time")
obj<-lapply(limits[2:length(limits)],function(x) lines(rep(x,2),
    par()$usr[3:4],lty="dashed"))

plot of chunk unnamed-chunk-4

or:

tree<-pbtree(n=26,tip.label=LETTERS)
tree<-paintSubTree(tree,29,"b","a",stem=0.5)
tree<-paintSubTree(tree,35,"c",stem=0.5)
colors<-setNames(c("black","blue","red"),letters[1:3])
plot(tree,colors)
xy<-markChanges(tree,plot=FALSE)
points(xy,pch=8)

plot of chunk unnamed-chunk-5

or:

## simulate data with multiple rates
x<-sim.rates(tree,setNames(c(1,10,20),letters[1:3]))
obj<-contMap(tree,x,plot=FALSE)
obj ## object of class "contMap"
## Object of class "contMap" containing:
## 
## (1) A phylogenetic tree with 26 tips and 25 internal nodes.
## 
## (2) A mapped continuous trait on the range (-12.257807, 1.528703).
plot(obj,sig=2)
xy<-markChanges(tree,plot=FALSE)
points(xy,pch=21,cex=1.5,bg="grey")

plot of chunk unnamed-chunk-6

obj$tree ## internal object of class "simmap"
## 
## Phylogenetic tree with 26 tips and 25 internal nodes.
## 
## Tip labels:
##  A, B, C, D, E, F, ...
## 
## The tree includes a mapped, 44-state discrete character with states:
##  198, 200, 218, 338, 366, 382, ...
## 
## Rooted; includes branch lengths.
plot(obj$tree,colors=obj$cols,lwd=4)

plot of chunk unnamed-chunk-6

More on this later!

No comments:

Post a Comment