Compute Hierarchy of an HBM

Share:

Description

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

Usage

1

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

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.

Value

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

Author(s)

Yoli Shavit

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

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
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)