R/partition.loss.dendrogram.R

In salso: Search Algorithms and Loss Functions for Bayesian Clustering

partition.loss.dendrogram <- function(estimate, x, loss) {
  isPSM <- isPSM(x)
  height <- partition.loss(estimate, x, loss)
  labels <- unique(estimate)
  specific.heights <- sapply(labels, function(label) {
    subset <- estimate == label
    y <- if (isPSM) {
      x[subset, subset, drop = FALSE]
    } else {
      x[, subset, drop = FALSE]
    }
    partition.loss(estimate[subset], y, loss)
  })
  min <- min(specific.heights)
  range <- max(specific.heights) - min
  candidates <- lapply(labels, function(label) {
    scale.factor <- if (range > 0) {
      (specific.heights[label] - min) / range
    } else {
      1
    }
    new.height <- height * scale.factor
    class(label) <- "dendrogram"
    attr(label, "members") <- 1
    attr(label, "height") <- new.height
    attr(label, "label") <- as.character(label)
    attr(label, "leaf") <- TRUE
    label
  })
  state <- list(estimate, candidates)
  while (length(unique(state[[1]])) != 1) {
    state <- agglomerate.dendrogram(state, x, loss)
  }
  state[[2]][[unique(state[[1]])]]
}

agglomerate.dendrogram <- function(state, x, loss) {
  estimate <- state[[1]]
  candidates <- state[[2]]
  tab <- t(apply(combn(unique(estimate), 2), 2, function(pair) {
    estimate[estimate == pair[2]] <- pair[1]
    c(pair, partition.loss(estimate, x, loss))
  }))
  which <- which.min(tab[, 3])
  pair <- tab[which, 1:2]
  estimate[estimate == pair[2]] <- pair[1]
  height <- tab[which, 3]
  candidates[[pair[1]]] <- {
    branch1 <- candidates[[pair[1]]]
    branch2 <- candidates[[pair[2]]]
    xx <- list()
    xx[[1]] <- branch1
    xx[[2]] <- branch2
    class(xx) <- "dendrogram"
    attr(xx, "members") <- attr(branch1, "members") + attr(branch2, "members")
    mid1 <- if (isTRUE(attr(branch1, "leaf"))) 0 else attr(branch1, "midpoint")
    mid2 <- if (isTRUE(attr(branch2, "leaf"))) 0 else attr(branch2, "midpoint")
    attr(xx, "midpoint") <- (mid1 + attr(branch1, "members") + mid2) / 2
    attr(xx, "height") <- height
    xx
  }
  candidates[pair[2]] <- list(NULL)
  list(estimate, candidates)
}