hbm (version 1.0)

hierarchy: Compute Hierarchy of an HBM

Description

hierarchy computes the hierarchy of a hierarchical block matrix computed with hbm.

Usage

hierarchy(hm)

Arguments

hm

a hierarchical block matrix computed with hbm, with:$hm_{i,j} =$ the minimal scale (iteration) at which i and j were clustered together, or 0 if i=j

Value

hierarchy returns a numeric value giving the hierarchy of the matrix.

Details

In a hierarchical matrix, computed with hbm, the behavior around the diagonal reflects the hierarchy of the association matrix. Specifically, for a hierarchical fractal-like structure we expect a non-decreasing series in the upper triangle of the matrix and a non-increasing series in the lower triangle. hierarchy counts the number of deviations from this behavior for each node: number of negative successive differences up to the diagonal and number of positive successive changes after the diagonal, and returns the negation of the mean number of changes across nodes.

References

hbm's website: http://www.cl.cam.ac.uk/~ys388/hbm/

See Also

generate.random.conf to see how to generate interesting chains hbm learn how to build hierarchical block matrices hbm's tutorials at http://www.cl.cam.ac.uk/~ys388/hbm/

Examples

set.seed(2)

n = 100 # chain size
#generate configurations 
conf = generate.random.conf(n, sd = 0.5, scale = FALSE)
#perturb the chain
conf.perturb.all = generate.random.conf(n, perturb = 1:n, sd = 0.5, scale = FALSE)
# and again with less perturbration 
conf.perturb = generate.random.conf(n, perturb = 10:50, sd = 0.5, scale = FALSE)

# compute the HBMs
hm.control = hbm(exp(-1*as.matrix(dist(conf))), 2)$hm
hm.perturb.all = hbm(exp(-1*as.matrix(dist(conf.perturb.all))), 2)$hm
hm.perturb = hbm(exp(-1*as.matrix(dist(conf.perturb))), 2)$hm

h.control = hierarchy(hm.control)
h.perturb = hierarchy(hm.perturb)
h.perturb.all = hierarchy(hm.perturb.all)
h = c(h.control, h.perturb, h.perturb.all)

# plot
plot(1:3, h, pch = 19, cex = 2, axes = FALSE, ylab = "Chain Hierarchy", xlab = "Condition")
axis(1, at = 1:3, labels = c("Control", "Perturbed-Partial", "Perturbed-All"))
axis(2)