Permalink
Browse files

new function & object class ctt

  • Loading branch information...
liamrevell committed Dec 4, 2017
1 parent b65922b commit f82043421c36fd216b6922f92c4109ddc6713af2
Showing with 199 additions and 4 deletions.
  1. +4 −4 DESCRIPTION
  2. +4 −0 NAMESPACE
  3. +157 −0 R/ctt.R
  4. +34 −0 man/ctt.Rd
View
@@ -1,6 +1,6 @@
Package: phytools
Version: 0.6-48
Date: 2017-11-29
Version: 0.6-49
Date: 2017-12-4
Title: Phylogenetic Tools for Comparative Biology (and Other Things)
Author: Liam J. Revell
Maintainer: Liam J. Revell <liam.revell@umb.edu>
@@ -13,6 +13,6 @@ ZipData: no
Description: A wide range of functions for phylogenetic analysis. Functionality is concentrated in phylogenetic comparative biology, but also includes a diverse array of methods for visualizing, manipulating, reading or writing, and even inferring phylogenetic trees and data. Included among the functions in phylogenetic comparative biology are various for ancestral state reconstruction, model-fitting, simulation of phylogenies and data, and multivariate analysis. There are a broad range of plotting methods for phylogenies and comparative data which include, but are not restricted to, methods for mapping trait evolution on trees, for projecting trees into phenotypic space or a geographic map, and for visualizing correlated speciation between trees. Finally, there are a number of functions for reading, writing, analyzing, inferring, simulating, and manipulating phylogenetic trees and comparative data not covered by other packages. For instance, there are functions for randomly or non-randomly attaching species or clades to a phylogeny, for estimating supertrees or consensus phylogenies from a set, for simulating trees and phylogenetic data under a range of models, and for a wide variety of other manipulations and analyses that phylogenetic biologists might find useful in their research.
License: GPL (>= 2)
URL: http://github.com/liamrevell/phytools
Packaged: 2017-11-29 12:00:00 EST
Packaged: 2017-12-4 12:00:00 EST
Repository:
Date/Publication: 2017-11-29 12:00:00 EST
Date/Publication: 2017-12-4 12:00:00 EST
View
@@ -134,6 +134,10 @@ S3method(print, phyl.pairedttest)
S3method(print, multi.mantel)
S3method(residuals, multi.mantel)
S3method(fitted, multi.mantel)
S3method(print, ctt)
S3method(print, multiCtt)
S3method(plot, ctt)
S3method(plot, multiCtt)
importFrom(animation, ani.options, ani.record, ani.replay, saveVideo)
importFrom(ape, .PlotPhyloEnv, .uncompressTipLabel, ace, all.equal.phylo, as.DNAbin, as.phylo, bind.tree, branching.times, collapse.singles)
View
157 R/ctt.R
@@ -0,0 +1,157 @@
## computing the mean number of character changes through time from a set of stochastic map trees
## written by Liam J. Revell 2017
ctt<-function(trees,segments=20,...){
if(!(inherits(trees,"multiSimmap")))
stop("trees should be an object of class \"multiSimmap\".")
tree<-as.phylo(trees[[1]])
changes<-sapply(trees,getChanges)
h<-max(nodeHeights(tree))
b<-segments
segs<-cbind(seq(0,h-h/b,h/b),
seq(1/b*h,h,h/b))
nchanges<-rep(0,b)
for(i in 1:length(changes)){
for(j in 1:length(changes[[i]])){
ind<-which((changes[[i]][j]>segs[,1])+
(changes[[i]][j]<=segs[,2])==2)
nchanges[ind]<-nchanges[ind]+1/length(changes)
}
}
LTT<-ltt(tree,plot=FALSE)
LTT<-cbind(LTT$ltt[2:(length(LTT$ltt)-1)],
LTT$times[2:(length(LTT$ltt)-1)],
LTT$times[3:length(LTT$ltt)])
ii<-1
edge.length<-rep(0,b)
for(i in 1:nrow(segs)){
done.seg<-FALSE
while(LTT[ii,2]<=segs[i,2]&&done.seg==FALSE){
edge.length[i]<-edge.length[i]+
LTT[ii,1]*(min(segs[i,2],LTT[ii,3])-
max(segs[i,1],LTT[ii,2]))
if(LTT[ii,3]>=segs[i,2]) done.seg<-TRUE
if(LTT[ii,3]<=segs[i,2]) ii<-if(ii<nrow(LTT)) ii+1 else ii
}
}
object<-list(segments=segs,nchanges=nchanges,edge.length=edge.length,tree=tree)
class(object)<-"ctt"
object
}
plot.ctt<-function(x,...){
h<-max(nodeHeights(x$tree))
args<-list(...)
if(!is.null(args$type)){
type<-args$type
args$type<-NULL
} else type<-"rate"
if(!is.null(args$show.tree)){
show.tree<-args$show.tree
args$show.tree<-NULL
} else show.tree<-FALSE
if(!is.null(args$add)){
add<-args$add
args$add<-NULL
} else add<-FALSE
if(is.null(args$ylim))
args$ylim<-if(type=="number")c(0,max(x$nchanges)) else
c(0,max(x$nchanges/x$edge.length))
if(is.null(args$xlim))
args$xlim<-c(max(x$segments),min(x$segments))
if(is.null(args$lwd)) args$lwd<-2
if(is.null(args$xlab)) args$xlab<-"time since the present"
if(is.null(args$ylab))
args$ylab<-if(type=="number") "mean number of changes" else
"mean number of changes / total edge length"
args$type<-"l"
args$x<-h-as.vector(t(x$segments))
args$y<-if(type=="number") rbind(x$nchanges,x$nchanges) else
rbind(x$nchanges/x$edge.length,x$nchanges/x$edge.length)
if(!add) do.call(plot,args)
else do.call(lines,args)
if(show.tree) plotTree(x$tree,add=TRUE,ftype="off",lwd=1,
color=make.transparent("blue",0.1),mar=par()$mar,
direction="leftwards",xlim=xlim)
}
sim.ctt<-function(tree,Q,anc=NULL,nmaps=100,...){
x<-as.factor(sim.history(tree,Q,anc=anc)$states)
flush.console()
cat("Starting stochastic mapping with simulated data vector.... ")
flush.console()
trees<-make.simmap(tree,x,Q=Q,nsim=nmaps,message=FALSE)
cat("Done.\n")
flush.console()
ctt(trees)
}
sim.multiCtt<-function(tree,Q,anc=NULL,nmaps=100,nsim=100,...){
object<-replicate(nsim,sim.ctt(tree,Q,anc=anc,nmaps=nmaps,...),simplify=FALSE)
class(object)<-"multiCtt"
object
}
getChanges<-function(tree){
states<-sort(unique(getStates(tree)))
nc<-sapply(tree$maps,length)-1
b<-which(nc>0)
nc<-nc[b]
xx<-vector()
H<-nodeHeights(tree)
for(i in 1:length(b)){
for(j in 1:nc[i]){
ss<-names(tree$maps[[b[i]]])[j+1]
x<-rep(H[b[i],1]+cumsum(tree$maps[[b[i]]])[j],2)
xx<-c(xx,setNames(x[1],
paste(names(tree$maps[[b[i]]])[j:(j+1)],
collapse="->")))
}
}
xx
}
print.ctt<-function(x,...){
cat("Object of class \"ctt\" consisting of:\n")
cat(" (1) a matrix (segments) with the beginning & ending time of each segment.\n")
cat(" (2) a vector (nchanges) with the mean number of changes in each segment.\n")
cat(" (3) a vector (edge.length) containing the total edge length of each segement.\n")
cat(" (4) an object of class \"phylo\".\n\n")
}
print.multiCtt<-function(x,...){
cat(paste(length(x),"objects of class \"ctt\" in a list.\n\n"))
}
plot.multiCtt<-function(x,...){
if(hasArg(alpha)) alpha<-list(...)$alpha
else alpha<-0.05
segments<-x[[1]]$segments
nchanges<-sapply(x,function(x) x$nchanges)
if(hasArg(type)) type<-list(...)$type
else type<-"rate"
if(type=="rate"){
edge.length<-sapply(x,function(x) x$edge.length)
}
obj<-list(segments=segments,nchanges=rowMeans(nchanges),
edge.length=rowMeans(edge.length),tree=x[[1]]$tree)
class(obj)<-"ctt"
lower<-floor(alpha/2*length(x))
upper<-ceiling((1-alpha/2)*length(x))
xx<-max(nodeHeights(x[[1]]$tree))-as.vector(t(segments))
xx<-c(xx,xx[length(xx):1])
y.lower<-if(type=="number") apply(nchanges,1,sort)[lower,] else
if(type=="rate") apply(nchanges/edge.length,1,sort)[lower,]
y.upper<-if(type=="number") apply(nchanges,1,sort)[upper,] else
if(type=="rate") apply(nchanges/edge.length,1,sort)[upper,]
y.lower<-as.vector(rbind(y.lower,y.lower))
y.upper<-as.vector(rbind(y.upper,y.upper))
yy<-c(y.lower,y.upper[length(y.upper):1])
args<-list(...)
if(!is.null(args$alpha)) args$alpha<-NULL
if(is.null(args$col)) args$col<-"blue"
if(is.null(args$ylim)) args$ylim<-range(yy)
args$x<-obj
do.call(plot,args)
polygon(xx,yy,col=make.transparent("grey",0.4),border=0)
}
View
@@ -0,0 +1,34 @@
\name{ctt}
\alias{ctt}
\alias{sim.ctt}
\title{Generates (or simulates) a 'changes through time' plot from a set of stochastic map character histories}
\usage{
ctt(trees, segments=20, ...)
sim.ctt(tree, Q, anc=NULL, nmaps=100, ...)
}
\arguments{
\item{trees}{an object of class \code{"multiSimmap"}.}
\item{segments}{number of segments to break up the history of the tree.}
\item{tree}{for \code{sim.ctt}, an object of class \code{"phylo"}.}
\item{Q}{for \code{sim.ctt}, a transition matrix to use for simulation.}
\item{anc}{ancestral state at the root node for simulation.}
\item{nmaps}{number of stochastic maps per simulation.}
\item{...}{optional arguments.}
}
\description{
This function generates a 'changes through time' plot in the style of a lineage-through-time (LTT) plot. It shows the mean rate or the mean number of changes per unit time from a set of stochastic character map trees.
}
\value{
An object of class \code{"ctt"} or \code{"multiCtt"}.
}
\references{
Revell, L. J. (2012) phytools: An R package for phylogenetic comparative biology (and other things). \emph{Methods Ecol. Evol.}, \bold{3}, 217-223.
}
\author{Liam Revell \email{liam.revell@umb.edu}}
\seealso{
\code{\link{ltt}}
}
\keyword{phylogenetics}
\keyword{plotting}
\keyword{comparative method}
\keyword{discrete character}

0 comments on commit f820434

Please sign in to comment.