###################<-60 Characters Wide->################### # part.recursive() # W. A. Green # Department of Geology, Yale University # walton.green@yale.edu # # Last Modified: May 6, 2005 # Function to partition a matrix of counts (into clusters or partitions) # that minimizes the log likelihood ratio of. # Depends on: part() #################### # Subroutines mat.overplot <- function(x, clus){ # Reorder matrix based on clus clusters <- unique(clus) clus.num <- length(clusters) x <- x[order(clus),] # Print out reordered matrix for(n in 1:nrow(x)){ for(m in 1:ncol(x)){ text(x = m/(ncol(x) + 1), y = 1 - (n/(nrow(x) + 1)), labels = paste(x[n,m]), col = sort(clus)[n]) for(o in 1:clus.num){ rect(0.05, 0.05, 0.95, 0.95) #rect(n/nrow(a) - ?, m/ncol(a) - ?, n/nrow(a) + ?, m/ncol(a) + ?) } } } if(any(clusters > 8)){ warning('colours may be repeated') } } # End of mat.overplot subroutine ######################################## # Body follows part.recursive <- function(x, col.struct = NULL, verbose = TRUE, ...){ if(is.null(col.struct)){ cat('no structure for the columns has been provided; clustering and reordering them automatically\n') col.clust <- hclust(dist(t(x))) col.struct <- merge2matrix(t(x), col.clust$merge, col.clust$order) } for(node in 2:nrow(col.struct)){ if(rowSums(col.struct)[node] < 1){ next }else{ x <- subx } # Perform the partitioning x.part <- part(x, verbose = 0, draw = 1) if(verbose){ print(x.part) } # ...and collect the clusters to reprocess clusters.to.reprocess <- names(table(x.part$cluster))[table(x.part$cluster) > 2] if(!(any(table(x.part$cluster) > 2))){ if(verbose){ cat('there are no clusters large enough to reprocess\n') } return(0) }else{ cat(length(clusters.to.reprocess), 'clusters with more than two rows:', clusters.to.reprocess, '\n') } submatrices.to.reprocess <- vector(length = length(clusters.to.reprocess), mode = 'list') for(i in 1:length(clusters.to.reprocess)){ submatrices.to.reprocess[[i]] <- x[x.part$cluster == clusters.to.reprocess[i],] } print(submatrices.to.reprocess) x.reordered <- x[order(x.part$cluster),] print(x.reordered) for(i in 1:length(submatrices.to.reprocess)){ index <- col.struct[node + 1,] print(node) print(index) subx <- submatrices.to.reprocess[[i]][,na.omit(index)] print(subx) quartz() subx.part <- part(subx, verbose = 0, draw = 0) } } return('hello') } # End of function