?methods
methods(plot)
?Reduce
help(rnorm)
require(graphics)#
#
dnorm(0) == 1/ sqrt(2*pi)#
dnorm(1) == exp(-1/2)/ sqrt(2*pi)#
dnorm(1) == 1/ sqrt(2*pi*exp(1))
par(mfrow=c(2,1))#
plot(function(x) dnorm(x, log=TRUE), -60, 50,#
     main = "log { Normal density }")
?rnorm
help.start()
x <- 1:30
x
require(ouch)
data(bimac)
data
bimac
rownames(bimac) <- bimac$node
tree <- with( bimac, ouchtree( node, ancestor, time/max(time), species  ) )
?ouchtree
plot(tree)
plot(tree, regimes = bimac[c("OU.3")])
plot(tree, regimes = bimac[c("OU.3")], lwd=3)
plot(tree, regimes = bimac[c("OU.LP")], lwd=3)
plot(tree, regimes = bimac[c("OU.1", "OU.3", "OU.4", "OU.LP")], lwd=3)
bimac[c("OU.3")]
bimac$OU.3
?brown
brown( log( bimac['size']), tree)
brown( log( bimac['size']), tree) -> bimac.bm
class(bimac.bm)
bimac.bm
h2 <- hansen( log(bimac['size']), tree, bimac['OU.1'], apha=1, sigma=1)
h2 <- hansen( log(bimac['size']), tree, bimac['OU.1'], alpha=1, sigma=1)
h2
plot( h2 )
h3 <- hansen( log(bimac['size']), tree, bimac['OU.3'], alpha=1, sigma=1)
h3
h4 <- hansen( log(bimac['size']), tree, bimac['OU.4'], alpha=1, sigma=1)
h5 <- hansen( log(bimac['size']), tree, bimac['OU.LP'], alpha=1, sigma=1)
h5
coef( h5 )
logLik(h5)
summary(h5)
unlist( summary(h5))
unlist( summary(h5)[c("aic", 'aic.c', 'sic', 'dof')])
brown( log( bimac['size']), tree) -> h1
h <- list( h1, h2, h3, h4 h5)
h <- list( h1, h2, h3, h4,  h5)
names(h) <- c("BM", "OU.1", "OU.3", "OU.4", "OU.LP")
h
sapply( h, function(x)  unlist( summary(h5)[c("aic", 'aic.c', 'sic', 'dof')]) )
print( sapply( h, function(x)  unlist( summary(h5)[c("aic", 'aic.c', 'sic', 'dof')]) ), digits=3)
print( sapply( h, function(x)  unlist( summary(h5)[c("aic", 'aic.c', 'sic', 'dof')]) ), digits=4)
print( sapply( h, function(x)  unlist( summary(x)[c("aic", 'aic.c', 'sic', 'dof')]) ), digits=4)
sapply( h, function(x)  unlist( summary(h5)[c("aic", 'aic.c', 'sic', 'dof')]) )  -> h.ic
h.ic
print( h.ic, digits=4)
sapply( h, function(x)  unlist( summary(x)[c("aic", 'aic.c', 'sic', 'dof')]) )  -> h.ic
print( h.ic, digits=4)
h5.sim <-  simulate( object = h5, nsim=10)
h5.sim
bimac
h5.sim
bootstrap( object= h5, nboot=10)
4.5+10.5+10.5
data(iris)
class(iris)
attributes(iris)
methods(class="data.frame")
methods(class="vector")
methods(class="numeric")
methods(class="matrix")
getAnywhere(plot.data.frame)
pairs
getAnywhere(pairs)
vignette("phylobase")
require(phylobaes)
require(phylobase)
data(geospiza)
class ? phylo4d
showMethods(class="phylo4d")
?plot(geospiza)
?pic
cat("((((Homo:0.21,Pongo:0.21):0.28,",#
   "Macaca:0.49):0.13,Ateles:0.62):0.38,Galago:1.00);",#
   file = "ex.tre", sep = "\n")#
tree.primates <- read.tree("ex.tre")
tree.primates
as(tree.primates, "phylo4")
as(tree.primates, "phylo4") -> phylo4tree
as(tree.primates, "phylo4") -> tree4
morph <- data.frame(row.names=labels(tree4), mass=rnorn( nTips(tree4)), mean=30, sd=10)
morph <- data.frame(row.names=labels(tree4), mass=rnorm( nTips(tree4)), mean=30, sd=10)
morph
morph <- data.frame(row.names=labels(tree4), mass=rnorm( nTips(tree4), mean=30, sd=10))
morph
tree4d <- phylo4d( tree4, tip.data=morph)
tree4d
plot(tree4d)
title("Phylo4d plot with default options")
plot(tree4d, center=F, scale=F, show.node.label=F, grid=F, ratio.tree=2/3, box=f)
plot(tree4d, center=F, scale=F, show.node.label=F, grid=F, ratio.tree=2/3, box=F)
quartz()
plot(tree4d)
morph2 <- data.frame(row.names=labels(tree4), mass=rnorm( nTips(tree4), mean=30, sd=10))
morph2 <- data.frame(row.names=labels(tree4), mass=rnorm( nTips(tree4), mean=30, sd=1))
morph2 <- data.frame(row.names=labels(tree4), mass=rnorm( nTips(tree4), mean=20, sd=10))
tree4d <- phylo4d( tree4, tip.data=c(morph, morph2)
)
morph
morph2
mm <- data.frame(morph, morph2)
mm
tree4d <- phylo4d( tree4, tip.data=mm)
plot(tree4d)
plot(tree4d, center=F, scale=F, show.node.label=F, grid=F, ratio.tree=2/3, box=F)
require(ouch)#
data(anolis.ssd)#
anolis.ssd
dat <- anolis.ssd[!is.na(anolis.ssd$species),]   # dropping the rows with NA (internal nodes)#
lssd <- dat$log.SSD         # copy the ssd data into a more convenient vector#
ecomorph <- dat$OU.7#
names(lssd) <- names(ecomorph) <- dat$species  # NOTE: names are required for matching in ape
ecomorph <- factor(ecomorph)
ssd.lm <- lm( lssd ~ ecomorph)#
anova( ssd.lm )
ssd.lm
setwd("/Users/marguerite/Docs_Work/Courses/Phylo4Dummies/chapters")
require(ape)#
atree <- read.nexus("Data/anolis.ssd.23tree.nex")  # read in as ape tree#
plot(atree)   # check that it looks OK
dat
ssdeco<- data.frame( lssd, ecomorph)#
ssd.pgls <- gls( lssd ~ ecomorph, correlation=corBrownian(phy=atree), data=ssdeco)#
summary(ssd.pgls)#
anova(ssd.pgls)   # produces the ANOVA F-ratio test for the overall effect of ecomorph
?lm
dat$lssd <- dat$log.SSD         # copy the ssd data into a more convenient vector#
dat$ecomorph <- dat$OU.7
ls
dat
dat$log.SSD
dat$log.SSD <- NULL
dat
dat
data(anolis.ssd)#
anolis.ssd
dat <- anolis.ssd[!is.na(anolis.ssd$species),]   # dropping the rows with NA (internal nodes)#
dat
dat[ c("species", "log.SSD", "OU.7" ]   # by column name#
dat[ c(2, 3, 7) ]   # by column number#
dat[ - c(1, 4, 5, 6) ]   # by excluding some column numbers (-)#
#
dat <- dat[ c(2, 3, 7) ]   # resave just the columns we need
names(dat) <- c( "species", "lssd", "ecomorph")   # rename all the columns#
names(dat)[2:3]  <- c("lssd", "ecomorph")     # rename just the 2, 3 columns (first was OK)
dat
lssd <- dat$lssd#
ecomorph <- dat$ecomorph#
names(lssd) <- names(ecomorph) <-  dat$species
plot( ecomorph, lssd)
ecomorph <- factor(ecomorph)
ssd.lm <- lm( lssd ~ ecomorph)#
ssd.lm#
anova( ssd.lm )
require(ape)#
atree <- read.nexus("Data/anolis.ssd.23tree.nex")  # read in as ape tree#
plot(atree)   # check that it looks OK
ssd.pgls <- gls( lssd ~ ecomorph, correlation=corBrownian(phy=atree), data=ssdeco)#
summary(ssd.pgls)#
anova(ssd.pgls)   # produces the ANOVA F-ratio test for the overall effect of ecomorph
?anolis.ssd
data(anolis.ssd)#
tree <- with(anolis.ssd,ouchtree(node,ancestor,time/max(time),species))#
plot(tree,node.names=TRUE)#
print(h1 <- brown(anolis.ssd['log.SSD'],tree))#
plot(h1)#
print(h2 <- hansen(anolis.ssd['log.SSD'],tree,anolis.ssd['OU.1'],alpha=1,sigma=1))#
plot(h2)#
print(h3 <- hansen(anolis.ssd['log.SSD'],tree,anolis.ssd['OU.7'],alpha=1,sigma=1))#
plot(h3)
print(h1 <- brown(anolis.ssd['log.SSD'],tree))#
plot(h1)#
print(h2 <- hansen(anolis.ssd['log.SSD'],tree,anolis.ssd['OU.1'],alpha=1,sigma=1))#
plot(h2)#
print(h3 <- hansen(anolis.ssd['log.SSD'],tree,anolis.ssd['OU.7'],alpha=1,sigma=1))#
plot(h3)
coef(h3)    # the coefficients of the fitted model#
logLik(h3)   # the log-likelihood value
summary(h3)   # everything in the print method except the tree+data
unlist(summary(h3)[c('aic', 'aic.c', 'sic', 'dof')])  # just the model fit statistics#
                                   # on a single line
h <- list(h1, h2, h3)   # store all our fitted models in a list
names(h)
h
h <- list(h1, h2, h3)   # store all our fitted models in a list
names(h) <- c("BM", "OU.1", "OU.7")
h
h[[1]]
sapply( h, function(x) unlist(summary(x)[c('aic', 'aic.c', 'sic', 'dof')]) )  #
                                   # table with all models
h.ic <- sapply( h, function(x) unlist(summary(x)[c('aic', 'aic.c', 'sic', 'dof')]) )  #
print( h.ic, digits = 3)
Sweave("IC_OUCHexamples.Rnw")
Sweave("IC_OUCHexamples.Rnw")
Sweave("IC_OUCHexamples.Rnw")
Sweave("IC_OUCHexamples.Rnw")
Sweave("IC_OUCHexamples.Rnw")
Sweave("IC_OUCHexamples.Rnw")
data(anolis.ssd)
anolis.ssd
?anolis.ssd
tree <- with(anolis.ssd,ouchtree(node,ancestor,time/max(time),species))
?ouchtree
plot(tree, node.names=T)
brown(anolis.ssd['log.SSD'],tree)
h1 <- brown(anolis.ssd['log.SSD'],tree)
h2 <- hansen(anolis.ssd['log.SSD'],tree,anolis.ssd['OU.1'],alpha=1,sigma=1)
anolis.ssd
h3 <- hansen(anolis.ssd['log.SSD'],tree,anolis.ssd['OU.7'],alpha=1,sigma=1)
plot( h3 )
coef(h3)
logLik(h3)
summary( h3)
h <-  list(h1, h2, h3)
names(h) <- c("BM", "OU.1", "OU.7")
sapply( h, function(x) unlist( summary(x)[c('aic', 'aic.c', "sic", "dof")]))
h.ic <- sapply( h, function(x) unlist( summary(x)[c('aic', 'aic.c', "sic", "dof")]))
print( h.ic, digits=3)