RiemStiefel (version 0.1.1)

st.hclust: Hierarchical Agglomerative Clustering on Stiefel Manifold

Description

Given the type of distance measure and agglomeration scheme method, gr.hclust performs hierarchical clustering on Grassmann manifold using fastcluster package, which returns the same object as stats package's implementation while providing more efficient computation. See hclust for more details.

Usage

st.hclust(
  x,
  type = c("intrinsic", "extrinsic"),
  method = c("single", "complete", "average", "mcquitty", "ward.D", "ward.D2",
    "centroid", "median"),
  members = NULL
)

Arguments

x

either an array of size \((p\times r\times N)\) or a list of length \(N\) whose elements are \((p\times r)\) matrix on Stiefel manifold.

type

type of distance measure; "intrinsic" or "extrinsic".

method

he agglomeration method to be used. This must be (an unambiguous abbreviation of) one of "single", "complete", "average", "mcquitty", "ward.D", "ward.D2", "centroid" or "median".

members

NULL or a vector whose length equals the number of observations. See hclust for details.

Value

an object of class hclust. See hclust for details.

Examples

# NOT RUN {
#-------------------------------------------------------------------
#      Generate a dataset with two types of Stiefel elements
#-------------------------------------------------------------------
#  group1 : first four columns of (8x8) identity matrix + noise
#  group2 : last  four columns of (8x8) identity matrix + noise

mydata = list()
sdval  = 0.05
diag8  = diag(8)
for (i in 1:10){
  mydata[[i]] = qr.Q(qr(diag8[,1:4] + matrix(rnorm(8*4,sd=sdval),ncol=4)))
}
for (i in 11:20){
  mydata[[i]] = qr.Q(qr(diag8[,5:8] + matrix(rnorm(8*4,sd=sdval),ncol=4)))
}

## try hierarchical clustering
#  compare 'intrinsic' and 'extrinsic' distance types
#  and use 'single' hclust option.
hint = st.hclust(mydata, type="intrinsic", method="single")
hext = st.hclust(mydata, type="extrinsic", method="single")

## visualize
opar = par(no.readonly=TRUE)
par(mfrow=c(1,2), pty="s")
plot(hint, main="intrinsic")
plot(hext, main="extrinsic")
par(opar)

# }