[R] ape comparative analysis query

From: Chris Knight <chris.knight_at_manchester.ac.uk>
Date: Thu 11 May 2006 - 00:16:16 EST


I've been comparing variables among objects (taxa) related by known trees, using phylogentically independent contrasts in the ape package, and want to move on to more complex models e.g. by using gls with appropriate correlation terms. My trees contain lots of (hard) polytomies and information about ancestors, which I've been includingcreating  fully dichotomous trees by using zero branch lengths. For instance:

library(ape)
tree<-read.tree(text="((B1:3,(((D1:5,C1:0):5,B2:0):4,(B3:3,B4:4):0):0):0,A1:0):0;")
plot(tree)

Where all the B taxa arise in a true polytomy from A1 and B2 is the ancestor of C1 which is the ancestor of D1

I have information on all the taxa and can create independent contrasts for this:

x<-c(B1=47,D1=43,C1=45,B2=50,B3=47,B4=48,A1=48)

y<-c(B1=2.9,D1=5.4,C1=2.8,B2=3.5,B3=3.2,B4=3.5,A1=1.8)
picx<-pic(x,phy=tree)
picy<-pic(y,phy=tree)

Which seems to make sense. However, if I try to use anything more sophisticated I run into problems:

>dat<-as.data.frame(cbind(x,y))
>model1<-gls(y~x,data=dat, correlation=corBrownian(phy=tree))
Error in corFactor.corStruct(object) : NA/NaN/Inf in foreign function call (arg 1)

or equally:

> model2<-compar.gee(y~x,data=dat, phy=tree)

[1] "Beginning Cgee S-function, @(#) geeformula.q 4.13 98/01/27"
[1] "running glm to get initial regression estimate"
[1] 15.6000 -0.2625
Error in gee(y ~ x, c(1, 1, 1, 1, 1, 1, 1), data = list(x = c(47, 43,  :
        NA/NaN/Inf in foreign function call (arg 14)

(NB using the 'data=' argument seems to be necessary in the gls case- if I don't, there is a further problem: 'Row names in dataframe do not match tree tip names. data taken to be in the same order as in tree. in: Initialize.corPhyl(X[[1]], ...)')

This seems to go away if I remove the root taxon (A1): tree2<-drop.tip(tree, "A1")
x2<-c(B1=47,D1=43,C1=45,B2=50,B3=47,B4=48) y2<-c(B1=2.9,D1=5.4,C1=2.8,B2=3.5,B3=3.2,B4=3.5) dat2<-as.data.frame(cbind(x2,y2))
model2<-gls(y2~x2,data=dat2, correlation=corBrownian(phy=tree2))

This raises various questions:

  1. Was I misleading myself that my independent contrasts were valid in the first place?
  2. What is it, if anything, about the root taxon that causes this issue, given that other taxa also have zero branch lengths?
  3. Is there any way of getting around this and including data on the root taxon, or am I better off just dropping it (ultimately I want to work with much larger trees (up to tens of thousands of taxa) where that one piece of information will become relatively less important)

Any help very much appreciated,

Chris

I'm working with ape 1.8-2 in R 2.1.1 under ubuntu 'Breezy' linux (unfortunately 2.1.1 is the latest easily available in breezy)

-- 
--------------------------------------------------------------
Dr Christopher Knight                      School of Chemistry
chris.knight@manchester.ac.uk     The University of Manchester
Tel: extn 64414                               Faraday Building
     +44 (0)161 3064414                              PO Box 88
Fax: +44 (0)161 3064556                       Sackville Street
                                            Manchester M60 1QD
 ` · . ,,><(((°>                                            UK

______________________________________________
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Received on Thu May 11 00:21:38 2006

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.1.8, at Thu 11 May 2006 - 14:10:06 EST.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.