Hi David,
collapse.singles() seems to work correctly if the "phylo" object is
correctly conformed. Some features of this class may seem annoying (and
Klaus and I alredy discussed about this), but these help for other
things. As a reminder the class "phylo" is defined here:
http://ape-package.ird.fr/misc/FormatTreeR_24Oct2012.pdf
Recently, I put a function on github to help code writters:
https://github.com/emmanuelparadis/checkValidPhylo
This could help you to detect problems that would be tough to find
otherwise.
Cheers,
Emmanuel
Le 12/06/2015 22:41, David Bapst a écrit :
Hi Klaus (and others),
Ah, I see! The real bug then appears to be in collapse.singles, as it
does not reorder the ID numbers in $edge. Here's a quick function to
return the root ID:
getRootID<-function(tree){
uniqueNode<-unique(tree$edge[,1])
whichRoot<-sapply(uniqueNode,function(x)
(sum(x==tree$edge[,2])==0))
rootID<-uniqueNode[whichRoot]
return(rootID)
}
And if we run it...
getRootID(tree) #original tree
[1] 151
getRootID(tree1) #after collapse.singles
[1] 151
getRootID(tree2) #after reorder.phylo
[1] 131
And we can see that although the number of nodes certainly changed
when collapse.singles was run, it didn't change the root node's ID.
So I guess this is really a collapse.singles bug report. Thanks for
the insight, Klaus!
-Dave
On Fri, Jun 12, 2015 at 2:27 PM, Klaus Schliep <[email protected]> wrote:
Hi David,
I found the bug. Somehow ape assumes on the one hand that the root is number
of tips + 1.
On the other hand the root of an phylo object is the node which is in
tree$edge[,1], but not in tree$edge[,2],
this is the case even in an unrooted tree. This annoys me for a long time.
The root of tree1 with the definition above is actually node 151 and not
131. Changing these solves your problem:
Here is a small function to do this for you:
switch.nodes<- function(tree, a, b){
ntips = length(tree$tip.label)
tree$edge[tree$edge==a] = 0L
tree$edge[tree$edge==b] = as.integer(a)
tree$edge[tree$edge==0L] = as.integer(b)
if(!is.null(tree$node.label)){
tmp<-tree$node.label[a-ntips]
tree$node.label[a-ntips] = tree$node.label[b-ntips]
tree$node.label[b-ntips] = tmp
}
tree
}
library(phangorn)
attr(tree1, "order") = NULL
getRoot(tree1)
tree2 = switch.nodes(tree1, 131, 151)
plot(tree2) # now works for me
Cheers and have a nice weekend,
Klaus
On Fri, Jun 12, 2015 at 2:53 PM, David Bapst <[email protected]> wrote:
Hello all,
As those of you who directly manipulate the guts of phylo objects in
your code (or construct new phylo objects whole cloth from
un-phylo-like data structures) have probably experienced, it is
sometimes easy to create $edge matrices that aren't accepted by ape
functions (I often use plot.phylo as my litmus for this).
When this occurs, various things can be done; I usually do the
following sequence:
collapse.singles()
reorder.phylo(,"cladewise")
read.tree(text=write.tree(,file=NULL))
...or something along those lines.
However, I've run into an issue where reorder.phylo() returns what is
essentially a scrap of the original $edge matrix, without warning.
Very worrisome! I'm using R version 3.2.0 and ape 3.3. Here's, my
script, including a call to a saved object on my website, so that the
error can be reproduced:
(Why do I have an .Rdata object listed with a .txt extension, you
might wonder? Because my school apparently sanitizes file extensions
on our hosted websites that it thinks are executables, archives or
doesn't recognize, but doesn't bother to check file innards when it
does recognize the extensions. Its easily hacked, at least.)
#################
library(ape)
load(url("http://webpages.sdsmt.edu/~dbapst/weirdTree_06-12-15.txt";))
#can I plot it?
plot(tree)
#nope
tree1<-collapse.singles(tree)
#any single-child nodes left?
sum(sapply(unique(tree1$edge[,1]),function(x) sum(x==tree1$edge[,1])==1))
#nope
#can I plot it?
plot(tree1)
#nope
#now reorder
tree2<-reorder.phylo(tree1,"cladewise")
tree2$edge
#now reorder with postorder
tree2<-reorder.phylo(tree1,"postorder")
tree2$edge
#################
As we can see, reorder.phylo with various ordering methods returns an
edge matrix with just six rows, from a tree that originally had
hundreds of edges. Most worrisome, it does this without any error
message. In this particular case, I retain the tree correctly if I
skip reorder.phylo and use the read.tree(write.tree) trick, but this
behavior of reorder.phylo is worrying.
Anyone have a clue what's going on here? If I am perhaps misusing
reorder.phylo(), is there an alternative approach for use in cleaning
phylo objects?
Cheers,
-Dave
--
David W. Bapst, PhD
Adjunct Asst. Professor, Geology and Geol. Eng.
South Dakota School of Mines and Technology
501 E. St. Joseph
Rapid City, SD 57701
http://webpages.sdsmt.edu/~dbapst/
http://cran.r-project.org/web/packages/paleotree/index.html
_______________________________________________
R-sig-phylo mailing list - [email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Searchable archive at
http://www.mail-archive.com/[email protected]/
--
Klaus Schliep
Postdoctoral Fellow
Revell Lab, University of Massachusetts Boston
_______________________________________________
R-sig-phylo mailing list - [email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Searchable archive at http://www.mail-archive.com/[email protected]/
