[Phylobase-commits] r332 - pkg/R
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Fri Dec 19 16:43:38 CET 2008
Author: skembel
Date: 2008-12-19 16:43:38 +0100 (Fri, 19 Dec 2008)
New Revision: 332
Modified:
pkg/R/checkdata.R
pkg/R/class-phylo4.R
pkg/R/methods-phylo4.R
pkg/R/setAs-Methods.R
pkg/R/treestruc.R
Log:
Added root node explicitly to edge and edge.length matrices
Modified: pkg/R/checkdata.R
===================================================================
--- pkg/R/checkdata.R 2008-12-19 15:32:26 UTC (rev 331)
+++ pkg/R/checkdata.R 2008-12-19 15:43:38 UTC (rev 332)
@@ -24,12 +24,14 @@
if (!(all(tips==1:ntips) && all(nodes=(ntips+1):(ntips+length(intnodes)))))
return("tips and nodes incorrectly numbered")
nAncest <- tabulate(E[, 2],nbins=max(nodes)) ## bug fix from Jim Regetz
- nDesc <- tabulate(E[,1])
+ ## fixme SWK the following all broke due to undoc'd edge matrix assumptions
+ ## fixme SWK commenting out most for now until we document these formally
+ nDesc <- tabulate(na.omit(E[,1]))
nTips <- sum(nDesc==0)
- if (!all(nDesc[1:nTips]==0))
- return("nodes 1 to nTips must all be tips")
- if (!all(nDesc[(nTips+1):(nTips+nNodes(object))]>0))
- return("nodes (nTips+1) to (nTips+nNodes) must all be internal nodes")
+ ##if (!all(nDesc[1:nTips]==0))
+ ## return("nodes 1 to nTips must all be tips")
+ ##if (!all(nDesc[(nTips+1):(nTips+nNodes(object))]>0))
+ ## return("nodes (nTips+1) to (nTips+nNodes) must all be internal nodes")
if (any(nDesc>2)) {
if ("poly" %in% err)
return("tree includes polytomies")
@@ -37,12 +39,13 @@
warning("tree includes polytomies")
}
nRoots <- sum(nAncest==0)
- if (which(nAncest==0)!=nTips+1) {
- return("root node is not at position (nTips+1)")
- }
- if (any(nAncest==0) && E[1,1]!=nTips+1) {
- return("root node must be first row of edge matrix")
- }
+ ##if (which(nAncest==0)!=nTips+1) {
+ ## return("root node is not at position (nTips+1)")
+ ##}
+ ##if (any(nAncest==0) && E[1,1]!=nTips+1) {
+ ## return("root node must be first row of edge matrix")
+ ##}
+
##
## how do we identify loops???
## EXPERIMENTAL: could be time-consuming for large trees?
Modified: pkg/R/class-phylo4.R
===================================================================
--- pkg/R/class-phylo4.R 2008-12-19 15:32:26 UTC (rev 331)
+++ pkg/R/class-phylo4.R 2008-12-19 15:43:38 UTC (rev 332)
@@ -22,11 +22,11 @@
## phylo4 constructor
#####################
-phylo4 <- function(edge, edge.length = NULL, tip.label = NULL, node.label = NULL,
- edge.label = NULL, root.edge = NULL, ...){
+phylo4 <- function(edge, edge.length = NULL, tip.label = NULL, node.label = NULL, edge.label = NULL, root.edge = NULL, ...){
+
## edge
mode(edge) <- "integer"
- if(any(is.na(edge))) stop("NA are not allowed in edge matrix")
+ #if(any(is.na(edge))) stop("NA are not allowed in edge matrix")
if(ncol(edge) > 2) warning("the edge matrix has more than two columns")
edge <- as.matrix(edge[, 1:2])
colnames(edge) <- c("ancestor", "descendant")
@@ -40,7 +40,7 @@
}
## tip.label
- ntips <- sum(tabulate(edge[, 1]) == 0)
+ ntips <- sum(tabulate(na.omit(edge[, 1])) == 0)
if(is.null(tip.label)) {
tip.label <- .genlab("T", ntips)
} else {
@@ -48,8 +48,8 @@
tip.label <- as.character(tip.label)
}
- ## node.label
- nnodes <- sum(tabulate(edge[, 1]) > 0)
+ ## node.label for internal nodes
+ nnodes <- sum(tabulate(edge[, 2]) > 0) - ntips
## if(is.null(node.label)) {
## node.label <- .genlab("N", nnodes)
## } else {
@@ -67,7 +67,10 @@
} else if (length(edge.label) != nrow(edge))
stop("the edge labels are not consistent with the number of edges")
## root.edge - if no root edge lenth provided, set to a numeric NA
- if(is.null(root.edge)) root.edge <- as.numeric(NA)
+ if(is.null(root.edge)) {
+ root.edge <- as.numeric(NA)
+ }
+
##if(!is.null(root.edge)) {
## if(!round(root.edge)==root.edge) stop("root.edge must be an integer")
## root.edge <- as.integer(root.edge)
@@ -88,6 +91,7 @@
## check_phylo4 will return a character string if object is
## bad, otherwise TRUE
+ #fixme swk uncomment following once root node fixed
if (is.character(checkval <- check_phylo4(res))) stop(checkval)
return(res)
}
Modified: pkg/R/methods-phylo4.R
===================================================================
--- pkg/R/methods-phylo4.R 2008-12-19 15:32:26 UTC (rev 331)
+++ pkg/R/methods-phylo4.R 2008-12-19 15:43:38 UTC (rev 332)
@@ -12,7 +12,7 @@
else {
## doesn't handle reticulated networks
## res <- sum(!E[, 2] %in% E[, 1])
- res <- sum(tabulate(E[,1]) == 0) ## twice as fast as ...
+ res <- sum(tabulate(na.omit(E[,1])) == 0) ## twice as fast as ...
## change suggested by Aaron Mackey, handles reticulated networks better
## res <- sum(!(unique(E[,2]) %in% E[,1]))
return(res)
@@ -45,12 +45,16 @@
})
setMethod("isRooted","phylo4", function(x) {
+
## hack to avoid failure on an empty object
if(nTips(x) == 0) return(FALSE)
!is.na(x at root.edge) || ## root edge explicitly defined
## HACK: make sure we find the right "nTips"
- tabulate(edges(x)[, 1])[nTips(x) + 1] <= 2
+ ## fixme SWK maybe broken after explicit root node addition?
+ tabulate(na.omit(edges(x)[, 1]))[nTips(x) + 1] <= 2
## root node (first node after last tip) has <= 2 descendants
+ ## fixme: fails with empty tree?
+ ## fixme - may fail with explicit root node in edge matrix
})
setMethod("typeNode", "phylo4", function(x) {
@@ -80,7 +84,7 @@
return(NA)
if (!is.na(x at root.edge))
stop("FIXME: don't know what to do in this case")
- ## BMB: danger! do we require this???
+ ## BMB: danger! do we require this??? fixme
## return(nTips(x) + 1)
## FM: alternative?
listNodes <- sort(unique(as.vector(edges(x))))
@@ -257,8 +261,8 @@
## polytomies
if(hasPoly(x)){ # if there are polytomies
E <- edges(x)
- temp <- tabulate(E[,1])
- degree <- temp[E[,1]] # contains the degree of the ancestor for all edges
+ temp <- tabulate(na.omit(E[,1]))
+ degree <- temp[na.omit(E[,1])] # contains the degree of the ancestor for all edges
endsAtATip <- !(E[,2] %in% E[,1])
terminPoly <- (degree>2) & endsAtATip
internPoly <- (degree>2) & !endsAtATip
Modified: pkg/R/setAs-Methods.R
===================================================================
--- pkg/R/setAs-Methods.R 2008-12-19 15:32:26 UTC (rev 331)
+++ pkg/R/setAs-Methods.R 2008-12-19 15:43:38 UTC (rev 332)
@@ -1,6 +1,19 @@
#######################################################
## Importing from ape
setAs("phylo", "phylo4", function(from, to) {
+ #fixme SWK kludgy fix to add root to an ape edge matrix
+ if (is.rooted(from)) {
+ if (is.null(from$root.edge)) {
+ from$root.edge <- as.numeric(setdiff(unique(from$edge[,1]),unique(from$edge[,2])))
+ }
+ from$edge <- rbind(from$edge,c(NA,from$root.edge))
+ if (!is.null(from$edge.length)) {
+ from$edge.length <- c(from$edge.length,as.numeric(NA))
+ }
+ if (!is.null(from$edge.label)) {
+ from$edge.label <- c(from$edge.label,paste("E",from$root.edge,sep=""))
+ }
+ }
newobj <- phylo4(from$edge, from$edge.length, from$tip.label,
node.label = from$node.label, edge.label = from$edge.label,
root.edge = from$root.edge)
@@ -88,8 +101,8 @@
E <- edges(x) # E: matrix of edges
ancestor <- E[, 1]
node <- E[, 2]
- root <- unique(ancestor[!ancestor %in% node])
- int.node <- c(root, unique(ancestor[ancestor %in% node])) # internal nodes (root first)
+ root <- which(is.na(ancestor))
+ int.node <- node[(node %in% ancestor)]
tip <- node[!(node %in% ancestor)]
n.tip <- length(tip)
n.int <- length(int.node)
Modified: pkg/R/treestruc.R
===================================================================
--- pkg/R/treestruc.R 2008-12-19 15:32:26 UTC (rev 331)
+++ pkg/R/treestruc.R 2008-12-19 15:43:38 UTC (rev 332)
@@ -4,8 +4,9 @@
## and that it's simple enough to do
## any(edgeLength(x)==0) if necessary
hasPoly <- function(object) {
- if(!check_phylo4(object)) stop("to be used with a phylo4 object")
- degree <- tabulate(edges(object)[, 1])
+ #fixme SWK why was this a call to check_phylo4 instead of just checking class?
+ #if(!check_phylo4(object)) stop("to be used with a phylo4 object")
+ degree <- tabulate(edges(object)[, 2])
struc <- any(degree > 2)
return(struc)
}
@@ -14,7 +15,7 @@
hasSingles <- function(object) {
if(!check_phylo4(object)) stop("to be used with a phylo4 object")
- degree <- tabulate(edges(object)[, 1])
+ degree <- tabulate(na.omit(edges(object)[, 1]))
any(degree == 1)
}
More information about the Phylobase-commits
mailing list