Re: [R-sig-phylo] 3. partial correlation with gls residuals? (Tom Schoenemann)
Hi Alejandro, I agree with almost all your points. A few additional comments below. Correlation gives the correlation among the estimated regression parameters: it is computed from the variance-covariance matrix extracted from the fitted object, model.fit$varBeta or vcov(model.fit) -- the latter is not to be confused with vcv(). It depends on the shape of the likelihood function and shows how much the estimation of a parameter may be affected by others. On phylogenetic signal, I think that AIC may be a suitable way to select the appropriate correlation structure: fit the same regression model with different corStruct and select the one with the smallest AIC. Note that often, the corStruct for a given regression model will not be the same than when considering each variable separately (e.g., gls(A ~ 1, ...) A more general comment to Tom: you seem to imply the following (causal) relationships (here '-' is not R's operator): C - A C - B It seems pretty clear to me that, if these relations hold, A and B will be correlated. But maybe you have other hypotheses/models in mind. Best, Emmanuel Alejandro Gonzalez wrote on 13/03/2012 16:04: Hello Tom, I will answer some of your questions below, Emmanuel and others may be better placed to reply to some of your questions: On 12, Mar 2012, at 7:49 PM, Tom Schoenemann wrote: Thanks Rob and Alejandro, OK, I did as suggested and ran a PGLS with A ~ B + C. I was hoping for some clarification of the actual results. Here is a summary: ** Generalized least squares fit by maximum likelihood Model: variableA ~ variableB + variableC Data: DF.B.A.C AIC BIC logLik -23.49499 -10.46914 16.7475 Correlation Structure: corPagel Formula: ~1 Parameter estimate(s): lambda -0.09862731 Coefficients: Value Std.Error t-value p-value (Intercept) 0.5229794 0.04740728 11.031625 0. variableB 0.2200980 0.05012508 4.390976 0. variableC 0.0620030 0.05128472 1.208996 0.2296 Correlation: (Intr) variableB variableB -0.813 variableC0.362 -0.837 Standardized residuals: Min Q1Med Q3Max -3.2951355 -0.4995948 0.2604608 0.8884995 2.8456205 Residual standard error: 0.2067425 Degrees of freedom: 100 total; 97 residual ** Are the following correct interpretations?: 1) Controlling for the phylogeny I used, variableB is associated with variableA independent of variableC, because the p-value of the beta weight for variableB is highly significant (0.) Yes, the model estimates a partial regression slope and standard error, that is the relationship between B and A controlling for the potential effects of C (that is potential relationship between C and A). 2) Controlling for the phylogeny I used, variableC is not significantly associated with variableA independent of variableB, because the p-value of its beta weight is not significant (0.2296) Yes, assuming there is no multicolinearity. See below... 3) There does not seem to be a significant effect of phylogeny on this relationship, since the ML lambda estimate is: -0.09862731 Yes, when lambda tends to 0, phylogenetic and non phylogenetic methods are expected to converge. I will note that an advantage of using PGLS methods is that the evolutionary parameter estimated simultaneously with model fit (in this case lambda) allows for adequate adjustment for phylogenetic signal. See a nice paper on this by Revell 2010 in Methods in Ecology and Evolution. Also, what exactly are the values listed in the Correlation section? Does the -0.837 entry indicate that variableB is correlated negatively with variableC controlling for variableA (and/or my phylogeny)? Maybe Emmanuel can answer this one. You could test whether B and C are strongly correlated, for example using PGLS. If they are repeating the analysis dropping B from your model should result in a significant correlation between C and A, indicating there is multicolinearity. Of course, the proper way to do all this is to examine your data prior to analyses. What are the two traits, maybe they are providing very similar information? Regarding my original plan to assess the independent relationship between variables using residuals, the Freckleton paper Alejandro kindly forwarded includes this comment: Note that to estimate the true slope for the effect of x2 using residual regression one would need to regress the residuals of the regression on y on x1 on the residuals of the regression of x2 on x1 (e.g. see Baltagi 1999, pp. 72ˆ74 for elaboration of this). (p. 544) This was what I remembered about the issue myself, though I haven't kept up on the literature Rob and Alejandro mentioned. However, Rob believes the residuals might not be independent of phylogeny, and
Re: [R-sig-phylo] Tree Permutation
Hi Sean,
Conceptually speaking it helps to think about tree permutations in terms of
snapping a tree at a random edge, and then re-attaching that edge somewhere
else on the tree. That saves the problems with picking two independent
nodes. There are standard routines for this kind of thing, particularly
nearest-neighbour interchange (NNI), subtree prune and regraft (SPR), and
tree bisection reconnection (TBR).
I know that NNI and SPR are implemented in Phangorn, but I couldn't find an
implementation of TBR. Others will know if there is one.
Cheers,
Rob
On 14 March 2012 23:55, Sean Downey s...@codexdata.com wrote:
Hi,
I've banged my head long enough with this one, so I hope someone can help.
I am trying to write a routine to permute the structure of a tree. Note,
this is different than phyloshuffle() in phylotools, which just shuffles
tips -- in contrast, I am trying to permute the structure of the tree. In
this case I cannot resample the original data because these trees are
manually constructed, yet I want a way or doing some statistics on them. If
there is a package doing this that I've missed? I would like to be able to
specify a number of iterations (easy) and the routine would repeatedly
switch random pairs of clades within the tree (not so easy, apparently). If
the routine is run for many iterations, the resulting tree would be
completely random. This seems simple, but given the dependencies in the
tree, you can't pick just any two nodes; for example, the second node
cannot be a descendant of the first node. Also, it seems like I need to
test that the two subtrees are structurally !
independent, but I'm not certain exactly what these criteria are. They
don't appear to have to be the same distance from the root in all cases.
The next problem is how to do this in ape. I'm reasonably familiar with the
phylo() object, but obviously since I'm posting this, I haven't been able
to make it work. Especially, the bind.tree() step.
Thanks much for any help,
Sean
## Tree Permutation
## Note, code does not work correctly!
library(ape)
library(distory)
atree-read.tree(text=((L1496,L1499),((L1498,L1490),((L1497,L1494),((L1493,L1492),(L1503,(L1502,(L1212,L1501)),L1489);)
atree-makeNodeLabel(atree, prefix = )
atree-compute.brlen(atree,1)
n=1 # No. itterations
plot.phylo(atree, show.node.label = TRUE)
ptree-atree # a copy to permute
for (i in 1:n){
# Pick the first node from the tree
a.node-as.numeric(sample(ptree$node,1))
a.node - 2 # testing
# This creates a subtree to sample the second branch from
extract.a-extract.clade(ptree,node=as.character(a.node))
drop.a.node-drop.tip(ptree,extract.a$tip.label, trim.internal=F) #trim F
to leave a place to attach other subtree
plot(extract.a, show.node.label = TRUE)
plot(drop.a.node, show.node.label = TRUE)
# Pick a node from somewhere ealse in the tree
b.node-as.numeric(sample(drop.a.node$node,1))
b.node - 5 #testing
# This creates a subtree to sample the second branch from
extract.b-extract.clade(drop.a.node,node=as.character(b.node))
drop.b.node-drop.tip(drop.a.node,extract.b$tip.label, trim.internal=F)
plot(extract.b, show.node.label = TRUE)
plot(drop.b.node, show.node.label = TRUE)
# drop.b.node contains whatever is left after clipping out the two
subtrees
# now I want to reassemble the tree switching extract.a and extract.b
bind.tree.1-bind.tree(drop.b.node,extract.a,where=b.node) # doesnt work?
plot(bind.tree.1, show.node.label = TRUE)
bind.tree.2-bind.tree(bind.tree.1,extract.b,where=a.node)
plot(bind.tree.2, show.node.label = TRUE)
#copy the results
ptree - bind.tree.2
}
#plot the permuted tree
plot(ptree, show.node.label = TRUE)
#END
*
Sean S. Downey, Ph.D.
Institute of Archaeology
University College London
31-34 Gordon Square
London
WC1H 0PY
s.dow...@ucl.ac.uk
http://www.homepages.ucl.ac.uk/~tcrnssd/
UK Cell: +44 (0)7542 41 0077
Belize Cell: +501 636-2939
___
R-sig-phylo mailing list
R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
--
Rob Lanfear
Research Fellow,
Ecology, Evolution, and Genetics,
Research School of Biology,
Australian National University
Tel: +61 2 6125 4321
www.robertlanfear.com
[[alternative HTML version deleted]]
___
R-sig-phylo mailing list
R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Re: [R-sig-phylo] Alternative to rapply() for nested lists of phylo objects?
Hi Dave,
I believe that this does what you want ('d'apply for Dave).
dapply - function(obj, FUN, class=phylo, default=NULL) {
ret - lapply(obj, function(x)
if (inherits(x, class)) FUN(x)
else if (is.atomic(x)) default
else dapply(x, FUN, class))
names(ret) - names(obj)
ret
}
library(ape)
tree-rtree(10)
## Named lists to check that things are where they should be.
list1-list(a=tree,b=tree,c=tree)
list2-list(e=list1,f=tree,g=list1)
dapply(list2, Ntip)
Basically this recursively applies maps itself over lists until it hits an
object of class 'class' (phylo by default). If it never finds an element of
that class by the time it reaches an atomic class (e.g., numeric) it will
return default (by default NULL):
list3-list(e=list1,f=tree,g=list1,h=NA)
dapply(list3,Ntip)
The above function won't be a good idea in all cases, and isn't much tested,
etc, etc. This also doesn't preserve class attributes, of the component lists,
but that could be done in the same way as names() is.
Cheers,
Rich
On 2012-03-14, at 2:39 PM, David Bapst wrote:
Hello all,
I am dealing with simulation results that produce nested lists of trees
(with unbalanced depths; i.e. the first list might have five lists of three
trees, the second element might have three lists of five trees, the third
element might just be a single tree). I want to do some analyses on each of
these trees, so one option would be flattening the list and recording their
original position in the list. Alternatively, I'd prefer to retain the
structure of the nested list. I recently discovered rapply() and thought
this was a perfect fix, only to discover that it appears to look at the
type() and not the class() of objects passed to it, so it ignores the phylo
objects as a whole and passes right on to their constituent elements.
For example:
library(ape)
tree-rtree(10)
list1-list(tree,tree,tree)
list2-list(list1,tree,list1)
rapply(list2,Ntip)
Playing with arguments in rapply doesn't stop this behavior; indeed, it
seems nothing will stop rapply from ignoring the type 'list' of the phylo
objects themselves. Does anyone know of an alternative to rapply that
considers the class of the nested objects instead?
Thanks,
-Dave B.
--
David Bapst
Dept of Geophysical Sciences
University of Chicago
5734 S. Ellis
Chicago, IL 60637
http://home.uchicago.edu/~dwbapst/
http://cran.r-project.org/web/packages/paleotree/index.html
http://home.uchicago.edu/%7Edwbapst/
[[alternative HTML version deleted]]
___
R-sig-phylo mailing list
R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
___
R-sig-phylo mailing list
R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Re: [R-sig-phylo] Testing for phylogenetic signal in proportions
It is really too few species to have a sufficiently powerful test for phylogenetic signal. See Blomberg et al. (2003), available on my website. Cheers, Ted Theodore Garland, Jr. Professor Department of Biology University of California, Riverside Riverside, CA 92521 Office Phone: (951) 827-3524 Wet Lab Phone: (951) 827-5724 Dry Lab Phone: (951) 827-4026 Home Phone: (951) 328-0820 Facsimile: (951) 827-4286 = Dept. office (not confidential) Email: tgarl...@ucr.edu http://www.biology.ucr.edu/people/faculty/Garland.html Experimental Evolution: Concepts, Methods, and Applications of Selection Experiments. 2009. Edited by Theodore Garland, Jr. and Michael R. Rose http://www.ucpress.edu/book.php?isbn=9780520261808 (PDFs of chapters are available from me or from the individual authors) From: r-sig-phylo-boun...@r-project.org [r-sig-phylo-boun...@r-project.org] on behalf of Rob Lanfear [rob.lanf...@gmail.com] Sent: Wednesday, March 14, 2012 8:18 PM To: r-sig-phylo Subject: [R-sig-phylo] Testing for phylogenetic signal in proportions Hi All, I've tried searching the literature for this but the search terms tend to give irrelevant papers, so I thought I'd ask for some guidance here. I have a set of 11 species, an ultrametric tree for the species, and for each species I have a set of proportional data for 5 chemical traits. E.g. for spp1 my data might be: chemicalA=10%, chemicalB=50%, chemicalC=1%, chemicalD=19%, chemicalE=20%. The set of traits for each species always sums to 100%, but I don't know, and can't measure, the absolute values of the traits. Biologically speaking, this is OK, because it's the proportions that I'm interested in. I want to test for phylogenetic signal in these traits, and estimate the rate of change of each proportion along the phylogeny. Can anyone point me to any appropriate references for the methods and pitfalls of attempting to do this with proportion data? I can see that proportional data will have some odd properties (non-independence of traits, bounded (i.e. non-Brownian) evolution, etc.), but the best way of accounting for these is not immediately apparent to me. Thanks, Rob -- Rob Lanfear Research Fellow, Ecology, Evolution, and Genetics, Research School of Biology, Australian National University Tel: +61 2 6125 4321 www.robertlanfear.com [[alternative HTML version deleted]] ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
