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
[[alternative HTML version deleted]]
_______________________________________________
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]/
