Hi, all. Last week I posted some code for extracting subtrees of trees in
hclust format. Petra Steiner quickly found an example for which the code
breaks, and sent it to me.
The problem seems to be that I had not considered the possibility of trees
with unlabeled nodes. In the new version of f.make.subtree (below), I steal
some code from plot.hclust to assign labels if there are none. That fixes
the problem. At the end I drop the labels again if the original tree didn't
have any.
Thanks to Petra for the help. If anyone has other examples, or other
suggestions for how to improve the code, please let me know.
Regards,
Matthew Wiener
Applied Computer Science & Mathematics Dept.
Merck Research Labs
Rahway, NJ 07090
732-594-5303
----------------------------------------------------------------------------
----
"f.make.subtree" <-
function (tree, groups, n = 1, add.element = NULL)
{
which.nodes <- which(groups == n)
names(which.nodes) <- NULL
## steal some code from plot.hclust
use.labels <-
if (is.null(tree$labels))
paste(1:(nrow(tree$merge) + 1))
else as.character(tree$labels)
subtree.labels <- use.labels[which.nodes]
subtree.order <- tree$order[use.labels[tree$order] %in%
subtree.labels]
which.merge.elements <- which(-tree$merge %in% which.nodes)
which.merge.rows <- sort(unique(which.merge.elements%%nrow(tree$merge)))
which.merge.rows[which.merge.rows == 0] <- nrow(tree$merge)
old.length <- 0
new.length <- length(which.merge.rows)
if (length(which.merge.elements) == 1) {
res <- list(merge = NULL, heights = tree$height[which.merge.rows],
order = 1,
labels = {if(is.null(tree$labels)) NULL else
subtree.labels},
method = tree$method,
call = tree$call, dist.method = tree$dist.method)
}
else {
while (new.length - old.length > 0) {
old.length <- new.length
new.rows <- numeric(0)
more.new.rows <- numeric(0)
row.elements <- as.vector(tree$merge[which.merge.rows,
])
pos.elements <- row.elements[row.elements > 0]
if (length(pos.elements) > 0)
new.rows <- pos.elements[apply(tree$merge[pos.elements,
, drop = F], 1,
function(x) {
any(x > 0) & all(-x[x
<
0] %in% which.nodes)
})]
more.new.rows <- which(apply(tree$merge, 1, function(x) {
all(x %in% which.merge.rows)
}))
which.merge.rows <- sort(unique(c(which.merge.rows,
new.rows, more.new.rows)))
new.length <- length(which.merge.rows)
}
merge.list <- tree$merge[which.merge.rows, , drop = F]
merge.height <- tree$height[which.merge.rows]
pos.elements <- sort(merge.list[merge.list > 0])
for (i in seq(along = pos.elements)) merge.list[merge.list ==
pos.elements[i]] <- which(pos.elements[i]
=which.merge.rows)
neg.elements <- merge.list[merge.list < 0]
minus.neg.elements <- -neg.elements
num.elements <- length(neg.elements)
change.table <- cbind(sort(minus.neg.elements), 1:num.elements)
for (i in 1:num.elements) {
merge.list[merge.list == -change.table[i, 1]] <- -change.table[i,
2]
}
old.order <- subtree.order
subtree.order <- match(use.labels[old.order], subtree.labels)
res <- list(merge = merge.list, height = merge.height,
order = subtree.order,
labels = {if(is.null(tree$labels)) NULL else
subtree.labels},
method = tree$method,
call = list(match.call(), tree$call), dist.method
tree$dist.method)
class(res) <- "hclust"
for (i in seq(along = add.element)) {
res[[add.element[i]]] <- tree[[add.element[i]]][which.nodes]
}
}
res
}
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-help-request at
stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
